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PREFACE 


The research described in this report forms the basis for 
a Ph.D. dissertation, "The Motion of Interconnected Flexible 
Bodies," by Arthur Stewart Hopkins, This dissertation serves 
as Volume 3 of the Final Report for Contract NAS 8-2 835 8, Mod, 6, 
sponsored by the NASA George C, Marshall Space Flight Center, 
and completed under the direction of Peter Likins, Principal 
Investigator. 
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ABSTRACT 


In this report the equations of motion for an arbitrarily inter- 
connected collection of substructures are derived. The substructures are 
elastic bodies which may be idealized as finite element assemblies. The 
substructures are subject to small deformations relative to a nominal 
state. Interconnections between the elastic substructures permit large 
relative translations and rotations between substructures, governed by 
Pfaffian constraints describing the connections. In the special case of 
screw connections (permitting rotation about and translation along a 
single axis) , constraint forces are eliminated and modal coupling is 
incorporated. 

This work is a generalization and extension of the recent literature 
directed at the. problem of flexible spacecraft simulation. It is an 
extension of Hurty's "component mode" approach in that interconnected 
elastic substructures are permitted large motions relative to each other 
and relative to Inertial space. It Is more general than the "hybrid 
coordinate" methods advanced by Likins and others, in permitting all sub- 
structures to be flexible (rather than only the terminal members of a 
topological tree of substructures). The presentation here is further 
distinguished from the spacecraft literature by its development of the 
basic relationships of continuum mechanics , on which spacecraft simula- 
tion models are based. 
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INTRODUCTION 


Prompted by increasing economic incentives for weight reduction 
and, in application to spacecraft, increasingly stringent pointing 
requirements of scientific instruments, there has been a growing inter- 
est in the dynamical behavior of flexible structures. Of course, the 
basic equations governing dynamical behavior were developed by scien- 
tists in the seventeenth and eighteenth centuries, principally by 
Newton, Euler and Lagrange. However, direct application of these 
equations to elastic continuum models is neither economically nor 
computationally feasible. Research in recent years has been directed 
toward obtaining efficient approximate techniques. Essentially two 
independent approaches were developed. 

The first approach was directed toward increasing computational 
efficiency for the detailed models required to accurately portray the 
small motion behavior of discretized models of flexible structures. 

This method is known as component modal analysis, and is generally ' 
attributed to Hurty (1965). As originally developed, the method 
applies to small motions of nonrotating structures. 

The second approach, in recognition of the need to include large 
relative motions between portions of the structure, was directed toward 
efficient computation of collections of point, connected rigid bodies. 
Early research in this area was conducted by Hooker and Margulies (1965) , 
followed by Roberson and Wittenberg (1966) using what may be classified 
as augmented body methods. A somewhat different approach to the same 
problem was developed by Velman (1967) and Russell (1969), who used 
what may be classified as nested body methods. As presented, both of 
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these methods apply only to point connected rigid bodies assembled in 
a configuration described as a topological tree. 

The failure of the rigid body methods to provide adequate defi- 
nition of flexible behavior, and the failure of the component modes 
methods to include large relative motions, led Likins (1967, 1968, 
1970) to a marriage of the two approaches. Using what he terms a 
hybrid coordinate approach, he has developed a procedure for treating 
the extremal bodies of a topological tree of point connected rigid 
bodies as flexible. As currently conceived, the method employs the 
Hooker-Margulies equations (as subsequently modified by Hooker (1970) 
to reflect the constraint force elimination techniques developed by 
Russell (1969)) for the central rigid bodies of the topological tree, 
and techniques similar to Hurty's for adding flexible extremal bodies. 

The present study consumates the union of the two approaches. In 
this study all of the bodies may be flexible, the configuration is 
arbitrary (topologically) and general differential (Pfaffian) con- 
straints are permitted. In addition to the general development, the 
special case of screw connections is treated and the resulting simpli- 
fications to the governing equations are detailed. 

The study is presented in three sections. The first section 
presents the theoretical' background and the derivation of the varia- 
tional equations governing the behavior of a material continuum in a 
rotating frame of reference. The second section applies the finite 
element method to obtain the governing ordinary differential equations 
for a substructure. These equations describe the small deformational 
motions of a structure idealized as a collection of finite elements. 
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This section also includes a demonstration of the adequacy of linear 
assumed displacement fields in the inherently quadratic Lagrangian for- 
mulation. A procedure for linearizing the equations about arbitrary 
frame motions is also presented, in conjunction with modal analysis 
procedures. The third section details procedures for connecting the 
substructures to form a structure, with the possibility of large rela- 
tive motions between substructures. In addition. to the equations for 
general Pfaffian constraints, the results are specialized to the case 
of screw connections (a rotation about and a translation along an axis) . 
In this special case (which includes most connections of practical use) , 
constraint forces are eliminated and the system equations are reduced 
by modal coupling. The summarizing fourth section presents the analysis 
procedure for an arbitrary dynamical system, as a means of abstracting 
the salient features of the study. 
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Section 1 


MATERIAL CONTINUA 

Newtonian mechanics provides an approximate description of a 
portion of the physical phenomena within man’s experience. Specifi- 
cally it accurately describes some of those phenomena for which 
relativistic effects are negligible; and whose spatial and temporal 
scales are large compared to atomic phenomena and small compared to 
cosmological phenomena. 

Physical phenomena are generally analyzed by means of the following 
procedure. First, a mathematical model is associated with the physical 
phenomenon being studied. Next, results are obtained from the model by 
purely mathematical methods. Finally, the mathematical deductions are 
extrapolated to the physical phenomenon. Experience has shown that 
these extrapolations accurately portray the physical phenomenon, if 
the model is developed with sufficient care. 

Since the second step in the procedure is mathematical, some 
understanding of mathematics is necessary. Mathematics consists of 
"mathematical objects" called "sets" (e.g. functions, numbers, geometri- 
cal objects) and "relations" between these objects (A.O). (To preserve 
the continuity of the presentation, the definitions of mathematical terms 
have been placed in Appendix A. When a term with a specific mathematical 
meaning is first used, it is set in quotation marks and followed by a 
reference to the appropriate section in the appendix.) The most funda- 
mental mathematical concepts cannot be defined and so are introduced, a 
priori, as primitive notions. The most fundamental rules governing 

JimCEDlNQ f AGE BLANKS W.TCnarf 
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mathematics cannot be deduced, and so are taken, a priori, as true 
"axioms" (A.O). 

Analogously, certain physical concepts will be accepted as primi- 
tive; and certain relationships will be taken as axiomatic. Physical 
concepts are, ultimately, high order abstraction of information obtained 
through the senses; therefore description of primitives is most conven- 
iently carried out in anthropomorphic terms. The purpose of this section 
is to identify the relevant physical concepts, and to associate a math- 
ematical model with them. Although a reasonably coherent development 
is intended, a mathematically rigorous treatment lies outside the 
scope of this study. 

1.1 Space and Time 

The concept of a physical system implicitly involves the notions 
space and time . Both will be taken as undefined primitives. For 
the range of phenomena being considered, space and time may be identi- 
fied with the usual intuitive concepts. Some of the more relevant 
characteristics of space and time will be listed. 

Space is infinite, homogeneous and isotropic. It is independent 
of time and physical phenomena. The concept of a point in space will 
also be primitive. A distance can be associated with any two points 
in space. An arbitrary reference distance (e.g. a meter) may be 
selected as a unit. The number of reference distances comprising the 
distance will be associated with the distance and the selected unit. 

Time is infinite, homogeneous and independent of space and 
physical phenomena. The concept of a point in time will be primitive. 
Time is anisotropic; for any two points, one is before and the other 
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after* A duration can be associated with any two points. An arbitrary 


reference duration (e.g. a second) may be selected as a unit. The 
number of reference durations comprising the duration will be associated 
with the duration and the selected unit. 

The determination of absolute position in space or time is pre- 
cluded by the infinite and homogeneous characteristics. Similarly, 
Isotropy precludes determination of absolute direction in space. 

An anthropomorphic primitive, the observer, will be introduced to allow 
tiie unique identification of points in space and time. The observer, 
as a spatially anisotropic collection of inhomogeneous points in space 
and time of finite extent, is able to uniquely identify the points of 
space and time by relating them to himself. 

The mathematical model associated with any specific, uniquely 
identified point in space will be a unique "mathematical object" (A.O), 
generically denoted by "s"; similarly, for any point in time, "t". The*' 
distinction between the physical and mathematical is customarily left 
implicit and the mathematical objects, s and t, are also called a point 
in space and a point in time. The model of space is the "set" (A.O) of 
all spatial points, which will be denoted "S",={s}, and called physical 
s P ace * The model of time is the set «f all temporal points, which will 
be denoted "T",- {t}, and called temporal space . It is convenient to 
introduce the "ordered pairs" (A.l) of spaLial and temporal points which 
will be called events and denoted "e",= (s,t). The set of all events, 
denoted "E", = {e}, is the "cartesian product" (A-l) of spatial and 
temporal points (E - S x T) and will be called event space . 
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The mathematical models associated with distance and duration will 
be "functions" (A.l), "a": S x S + R and "t": T x T -*■ R, assigning non- 
negative real numbers to pairs of points in physical space or temporal 
space. The functions, a and t, define "metrics" (A. 2) on S and T. The 
function, "e E x E + R, p >_ 1, defined by 

epKjS.jt), ( 2 s, 2 t)] = ([a(jS, 2 s)] P + [T( 1 t, 2 t)J P y L/p , (1.1-1) 

is a metric for E. Similarly, in the limit as p 03 , 

e co [( 1 s, 1 t) , ( 2 s, 2 t)] = max(a( 1 s, 2 s) , x( 1 t, 2 t)J (1.1-2) 

is also a metric. All of the £ -metrics are "strongly equivalent" (A. 2). 

P 

Physical space, S, with a, temporal space, T, with T, and event space, 

E, with any of the £^, are all "metric spaces" (A. 2). 

If the "open sets" (A. 3) of S,T or E are defined as arbitrary 

unions of "open balls" (A. 2) the collections of all open sets are 

"topologies" (A. 3), specifically they are "metric topologies" (A. 3); 

and the sets, along with the topologies, define "topological spaces" 

(A. 3). The topology on E obtained from any of the £ -metrics is the 

1 ? 

same as the "topological product" (A. 3) of S and T. 

1.2 Coordinates 

The mathematical models of space and time will be related to the 
real numbers to facilitate manipulation. The observer can establish 
this relationship by constructing a coordinate system. He may 
specialize an event as the origin, and certain spatial points as 
dextral, orthogonal, coordinate axes. The duration between an event 
and the origin, and the distances between the event's geometric 
projections on the axes and the origin, along with appropriate sign 
conventions, establish the usual coordinate system . 
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The mathematical model of a "coordinate system" (A. 4) is a 

function from the topological space onto a cartesian product of the 

3 4 

real numbers (e.g. S ■+• R , p: T -+• R, or ]i: E R ) . Taking this 

b i t 

function as the only "chart” (A. 4) in an "atlas" (A. 4), physical, 
temporal and event space become "C manifolds" (A. 4). The product 
topology, E, with the chart, " y"» = is called the "product 

manifold" (A. 4). 

The manifold provides a mathematical model that transcends the 

coordinate system. Since a manifold, by definition, includes all 

"admissible charts" (A. 4), a change in coordinates does not change the 

manifold. For physical space, origin translations and axes rotations 

do not change the manifold; polar, spherical and the other customary 

curvilinear coordinates (with minor restrictions) are all just special 

cases of the admissible coordinates. 

The metrics, a: S * s -+■ R and x: T x t -+■ R, have a simple form in 

3 12 3 

terms of the usual coordinates, R , ' = (x ,x ,x ), and 

4 

T y: T R, = (x ); namely 


a( l S ’2 s) = (t xl (i s ) " x 1 ( 2 s)] 2 + [x 2 ( 1 s)~ x 2 ( 2 s)J 2 + tx 3 ( lS )- x 3 ( 2 s)] 2 ) 

( 1 . 2 - 1 ) 

and 

x( 1 t, 2 t) = j x^t) - x^( 2 t)| . (1.2-2) 


1/2 


1.3 Hatter 

Certain events, called material events , are characterized by the 
presence of matter . Matter will be taken as an undefined primitive 
notion, some of the more important characteristics of which will be 
listed. 
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Matter occurs as aggregates of spatial and temporal points, never 
as an isolated point in space or time. Mathematically, material events 
occur in "open sets" (A. 3), Only aggregates of finite spatial and 
temporal extent will be treated, so the open sets will be "bounded" (A. 2) 
(and their closures will be compact). An aggregate of matter, over 
which the physical characteristics vary smoothly, will be called a 
material continuum and will be denoted "M". The characteristics are 
called material properties and are determined by the constituent atoms. 
Mathematically, the material properties are functions, and varying 
smoothly corresponds to being a "C^ Map" (A. 4). 

Inertia is an important characteristic of matter that will be 
taken as primitive. It can be measured by the material property mass 
density . An arbitrary reference mass density (e.g. a kilogram/meter 3 ) 
may be selected as a unit. The number of reference mass densities 
comprising the mass density will be associated with the mass density 
and the selected unit. Mathematically, the mass density is the C°° 
map "p": M -+■ R assigning a non-negative real number to each point in a 
material continuum. 

Another important characteristic of matter is that it interacts 
with other matter. This interaction is called force , and may be 
associated with the anthropomorphic notion of a push or pull. Force 
will be taken as a primitive, even though it is intimately related to 
inertia. Forces may generally be categorized as either long range 
(remote or field) forces or as short range (contact) forces; depending 
on whether the distance between the interacting matter is very large 
or very small in comparison to the dimension of the material continuum. ' 
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Long range forces can be measured by body force densities . An 


3 

arbitrary reference body force density (e.g. a Newton/meter ) may be 
selected as a unit. The number of reference body force densities com- 
prising the body force density will be associated with the body force 
density, and the selected unit. If the material properties are constant 
over the material continuum, the long range forces vary slowly over the 
continuum. Gravity is a primary example of a long range force. 

Conversely, short range forces vary extremely rapidly over the 
material continuum. If the interaction between matter exterior to and 
interior to a material continuum is considered, the forces are extremely 
high in the region of the "boundary" (A. 3) and drop to negligible values 
at an extremely short distance into the "interior" (A. 3), It is custom- 
ary to treat the integral of these forces, normal to the surface, as a 

surface force density defined on the boundary. An arbitrary reference 

2 

surface force density (e.g, a Newton/meter ) may be selected as a unit. 
The number of reference surface force densities comprising the surface 
force density will be associated with the surface force density and 
the selected unit. Forces at the points of contact between bodies are 
primary examples of short range forces. , 

Forces are characterized by direction as well as magnitude. The 
mathematical model for the surface or body force density at a point in 
a material continuum is a "vector" (A. 6) which is a member of a "vector 
space" (A. 6) at the point. The "tangent space" (A. 5) at a point in a 
manifold is a natural vector space to associate with a manifold. The 
forces at each event of a material continuum collectively form a 
"vector field" (A. 5). 
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The final characteristic of matter to be described here is 


traceability . Matter persists through time; exactly the same matter 
is present at any point in time in a material continuum (i.e., matter 
neither appears in nor disappears from a material continuum). For any 
specified event in a material continuum, and any (possibly different) 
given time, there is exactly one event in the material continuum which 
occurs at the given time and which may identified as being the identical 
matter (i.e,, consisting of the same atoms) as the previously specified 
event. 

The subset of the events in M for which the time "proj ection" (A.4) 
is some reference time "t " will be called the time slice at t . It is 

o 0 

a "submanifold" (A. 4) of M (which is an "open submanifold" (A. 4) in E) . 

A "diffeomorphism 11 (A, 4) may be selected from the time slice onto the 
subset "S" of physical space. The mapping is usually just the spatial 
projection, but it may be used to establish a reference spatial state 
of the material which does not occur in the time period being treated 
(e.g., an undeformed state). The image, "s, is an open submanifold of S, 
and will be referred to as material space . The temporal projection of 
M is an open submanifold of T, which will be denoted "T,” and may be 
referred to as the time range of the material continuum. The product 
manifold, "M", = S x T, will be called material event space . It may be 
noted that M and M are different and provide two different means of 
identifying matter. The first, M, by its reference position and the 
time (referred to as Lagrangian) and the second, M, by its actual 
position at the given time and the time (referred to as Eulerian) . 
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The property of traceability can now be formulated as the 
diff eomorphism, "ft": M -+ M, which assigns to each material point (i.e., 
the reference position of the matter), "s" £ S, and to each time, 

"t" € T, that is to each material event, "m" € M, the event in the 
material continuum, "m" € M, which represents the location at time t 
of the matter occupying s in the reference state. The event, m, will 
be called the position of the matter. It may be noted that this 
implicitly limits the material continua to those subsets of event space 
for which there is a 1-1 correspondence between the spatial points in 
any time slice and the spatial points in the reference state. 

The manner in which M is constructed implies certain relationships 

for II. The times are the same (e.g., if (s,t) ■- H(s, t) , then t » t). 

The diff eomorphism from the time slice to S is just the restriction to 

t = t of the inverse of II (e.g,, (s, t ) “ II ^ (s,t )). The usual 

r — i i — * 

coordinate system, lx j, on M is the product of the restrictions of y 

s 

and to S and T. In terms of the usual coordinates, the "coordinate 
expression" (A, 4) for ft is 

t rr X / 1 3 4 \ / i rj <1 V ’ 

xoH=Tr(x,x,x,x) . (1.3-1) 

1.4 Velocity 

For any specific s € S, the position function assigns a point in 
M to each time, t, 6 T. These points in M or M form a "differentiable 
curve" (A. 5) parameterized by the usual time coordinate. These curves 
represent the trajectory of the matter which occupies s in the refer- 
ence state. The "tangent vector" (A. 5) to these curves at any event 
has a physical interpretation as the velocity at the event. The 
collection of all such tangents is modeled as a "vector field" (A. 5) 
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called the velocity field and denoted "V" or "V". Each of the curves 


is the "integral curve" (A. 5) that passes through H (s, t Q ) or (s, t Q ) , 
respectively. 

— -1 —1 — — 2 — 

An integral curve in M has the form "y(u)" = y (x (s),x (s), 

x^(s),u); that is it is the fourth "coordinate curve" (A. 5). Its 

^ _ 

tangent is "y^m)", = — -r (m) . The velocity field, V, in M is thus 

— 9 

V = — t ♦ The tangent to the integral curve in M, II o y (u) , is 
9x 


obtained from the "differential" of II (A. 5) as 


V[H(m)] = n*[y* <m) ] = ^ -2 ^ (m) ^-roll(m) 

9x 9 k 1 


(1.4-1) 


Recalling the coordinate expression for H, (1.3-1), the velocity field 


in M is 


.. SCx^JI) „-l 3 Sir 1 - .-1 S 

V = — t'-oH ~T = — royoU t 

dx d.X 9u v dX 


(1.4-2) 


It may be noted that since x o II = x the component of — -r is 1. (In 

4 1 2 3X 3 

the usual notation with x = t, x oll = r,xol[=r,xoll = r, 

* x y Z* 


— r- - i , — — = j, - — - = k , this becomes 
9x 9x 9x J 

9r- 9r A 9r A „ 

V 5- i + — i f + — * k + . 

9t 9t j at at * 


(f.4-3) 


the spatial part of which is the usual definition of velocity.) 

The velocity at any point, m, is algebraically a vector in the 


vector space the "tangent space" (A. 5) to the material continuum. 

I 9 

The "coordinate vector fields," (A. 5) "j r ", form a "basis" (A. 6) 

1 ^x 1 

for the space. Among all the bases, the coordinate vector fields 


14 



satisfy the additional relationship of the "usual tensor transform" 


(A. 5). If "{y }" are also coordinates, then V = a r (m) - b 1 r (m) 

i 1 ax 1 3x1 3yl 

implies a * b J — — (m) . The velocity at a point is, thus, a "tensor” 

Sy^ 

of "type (1,0)", a "contravariant vector" (A. 7), belonging to the 


"tensor space" M (A. 7). 

m 

1.5 Reference Frames 


The notion of a reference frame is equivalent to the concept of 
the event space of an observer. If there are two observers, they will 
be, in general, in motion relative to one another. The spatial points 
that appear fixed to one observer appear to be moving relative to the 
other observer (i.e. , the observers' event spaces or frames of refer- 
ence are different). If at some point in time the two observers 
mutually agree upon the identification of points in space and time, and 
upon a coordinate system, their temporal spaces and coordinates will 
be the same. However, at any time other than the reference, their 
physical spaces may differ (e.g. a point in space that they agreed was 
the same at t Q , they will, in general, interpret as different points 
at' any other time). If they both select the reference time for the 
time slice and the same reference state, material event space will be 
the same for both observers. 


The mathematical model of this relationship is a dif feomorphism, 

A : E’ -> E, (1.5-1) 

which assigns to any event in the second observer’s space, the event in 
the first observer’s space that they would mutually agree is the 
identical event physically. Let the p’ - y coordinate expression for 
this relationship (1.5-1) be 
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i . .i, '1 '2 '3 ’4, 

x o A = A (x ,x , x , x ) 


(1.5-2) 


i 1 i 

for the usual coordinates, tx } and {x } , on E and E 1 . 


Since the distance and duration functions are associated with 
space and time, and are not dependent on the observer, they should 
appear the same in both observers' physical and temporal spaces. The 
metrics are equivalent to "bilinear forms" (A. 7) on the manifold. In 
fact, a and T may be obtained as the "distance" (A. 7) between points 
with respect to the bilinear forms (in terms of the usual coordinates) 


dx a ® dx^ 


(1.5-3) 


4 4 

dx ® dx 


(1.5-4) 


(The Greek postscripts indicate a range of 3, a, 3 - 1,2,3. Latin 
postscripts will be used when the range is 4.) These forms are 
"A-related" (A. 8) to the bilinear forms 

10^1 3(X S o H ') dx ,Y @dx ,S (1 

“ S 9x' Y 3x' 6 U 


dx ,4 @dx ,4 _ 

For the bilinear forms to have the expressions on E 1 , 
dx ,a ® dx'^ 


4 4 

dx’*®dx’* , 


the X must satisfy 


3(X“ o u’) 3(X a on’) r 

n _ . ~ 0 , 
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and 


3(X 4 o y') 3{X 4 o y’) 

f I 1 


,4 


.,4 


(1.5-10) 


3x' 

throughout E 1 . For these relationships to hold at every event in E', 
the Indicated partial derivatives must be constants, and consequently 
the associated must be linear in the appropriate coordinates. The 
p' - y coordinate expression must, therefore, have the form 

x^ o A = r\x' 4 ) + (x’ 4 ) x- 1 + c 1 ^ (x* 4 ) x'^ + c T ^ (x* 4 ) x’\ 

2 a 2 r . 4 x . ,2 t , 4 n ,1 x , 2 , , 4 . ,2 t 2 , t 4 N T 3 

x o A = r (x ) + c ^ (x ) x + (x* ) x' + c’^ (x ) x T , 

3 „ a 3, f4x . ,3 , ,4 S ,1 , 3 , ,4, ,2 t 3 / ,4, t 3 

x o A = r (x T ) + c ^ (x ) x + c ^ (x ) x + c r ^ (x’ ) x , 

x 4 o A = x' 4 (1.5-11) 

and satisfy the constraints 


• Ci / 1 4. i ot f ,4 n 

c’ (x f ) c’ (x' ) - 6 a . 

B y By 


(1.5-12) 


The arbitrary translation in time coordinates has been set to aero, so 
the observers share the same time coordinate. 


The coefficients of this expression have a physical interpretation. 
The r terms represent the unprimed components of the displacement from 
the unprimed to primed origin. The c'^ terms represent the direction 
cosines between the x and the x ,p axes. 

The transform is apparently determined by the twelve real valued 
functions of time, r and c 1 ; however, the constraints limit the number 
of independent c' to three. There are several alternative parameteri- 
zations of the c’ part of the transform. The Euler-Rodrigues parameters 
are an example of a four parameter set, and the Euler angles are. an 


17 



example of a three parameter set. All three parameter sets have 


singularities (i.e., transforms that cannot be represented by the set). 
The inverse map, A 1 : E *+ E', has a p - U* coordinate expression 

x 1 ^ o A * - - r T *(x 4 ) + c^(x 4 ) x 1 + c^x 4 ) x 2 + c^(x 4 ) x\ 

,2 “ A -l ,2 , 4. 2 . 4, 1 2, 4. 2 ( 2 . 4, 3 

x o A = - r (x ) + c^(x ) x + c 9 (x ) x + c^(x ) x , 

,3 .-1 , 3 , 4. 3,4 1 3.4. 2 , 3, 4. 3 

x o A - - r (x)+ c^(x ) x + c 2 ^ x ) x + c 3 ^ x ) x 


and 


,4 A -1 4 

x o A - x , 


subject to the constraints 

oc / 4 ■. a , 4. , 

c e < x > C Y < x > - 6 gy • 


(1.5-13) 


(1.5-14) 


-1 


Since the composition of A and A must be the identity map, 

a a, 4. , ,a . 4. , ,6. 4. , , 4 S y., 

x = r (x ) + c' fl (x ) [- r' (x ) + C (x ) x r ] 


Y 


= r a (x 4 ) - c*“ (x 4 ) r’^(x 4 ) + c‘?J (x 4 ) (x 4 ) x Y . (1.5-15) 


B ‘ v ~ 7 ’ w 3 ' ' Y 


From this it may be concluded that r, r 1 , c, and c* satisfy the relations 
a / 4. ,a ,4. ,g . 4. 

r ( x ) ” c’p (x ) r’ (x ) (1.5-16) 

and 

C ’e (x4) c y (x4) = 6 Y * (1.5-17) 

Recalling the constraints (1.5-12 and -14), it may be seen that 

(1.5-18) 

that is, that C is the transpose of C 1 . Since its transpose is its 
inverse, C is an orthogonal matrix. 

A contravariant vector field, V’ :M f TM T , (e.g. velocity) has a 
A-related vector field, V : M -»■ TM, defined by 


c a (x4) = c 'b (x4) : 
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V = A* (V * ) 


- 1 0 A “1 9(x J o A) A -1 




o A 


(1.5-19) 


9x’^ 9x j 

(for the usual coordinates) . Recalling the coordinate expression for 
A (1.5-11), this is 


V - |v’ 3 c’^Cx’ 4 ) + v’ 4 

p 


+ v 4 o A- 1 -A 


3r a (x’ 4 ) 3 c 'b( x ' 4 ) 


,.e 


3x 


,4 


9x 


,4 


o A 


-1 9 


3x“: 

(1.5-20) 


9x 




where the v* are the components of V T with respect to the usual coordi- 

* 4 4 

nates. Note that v 1 = v = 1 for velocity. - 

As an alternative procedure, the real valued functions on M, 

x' 1 o A may be taken as coordinates, iy^"}, on M. The expression for 

V under a change of coordinates is provided by the usual tensor 

;,i 9 4 ^ 


transformation V =v 


,i 3x J 9 
= v — 


9y‘ 


9y X 9x j 


So 


V = lv ,Q C^ (y 4 ) + V * 4 

a 7 


.3/4, 


3r^Q_ + ^_a a 

l 3y 4 3y 4 


.8 


9x" 


+ v 


,4 9 


(1.5-21) 


9x 


Thus reference frames may be considered either as special classes of 
diff eomorphic manifolds, or as special classes of coordinatizations of 
a manifold. 

It remains to be shown that the transform of the primed observer's 
velocity is actually the velocity that the unprimed observer would 
assign. If 11' : M ->■ M* is the mapping from material event space to the 
primed material continuum, then II = A o H' , since M is the same for both 
observers. If the coordinate expression for IT' is 
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(1.5-22) 


,i , ,i 1 ~2 —3 —4. 

x o K = tt t (x , x , x , x ) , 

then the velocity in M’ is 


V - l^^L. on 1 ' 1 - — T - “'I 1 
3x Sx* 1 3u 


97T — nt -l 9 
o y o IT — 


3x 


,i * 


(1. 5-22b) 


Substituting II = A o H' in x^ o II = tt'*' o p yields x^oAoTI’ = o y. 
Inserting the coordinate expression for A (1.5-11) yields 


iT a o y = (r a (x f4 ) + c’^x* 4 ) x^] o II 1 and 

4 — . ,4 rr t 

7T o yi - X O n . 

The usual velocity in M is 


(1.5-23) 


V = o if 1 i-r = m „ n-1 3— 

ai 4 ax 1 a^ 4 ax 1 

Substituting (1.5-23) in (1.5-24) yields 


+ c*' e o r, + c •“(!*> 

aP ai 4 e aJ 4 


(1.5-24) 


o n 


-1 a 


a 


3x 


3x 

(1.5-25) 


1 4 


—4 


since x’ " o II* - x - ’. Noting that x 4 o iT 1 « x 4 and IT o IT _i = A' 1 , 


r-1 .-1 


V = ( [c^(x 4 )][v t3 o A"" 1 ] + 4- 

3x 

3 


rv | Ol / 4. 

3c g(x ) 


1 3x 


[x’^oA- 1 ! U- 


3x“ 

(1.5-26) 


as was obtained before (1.5-20 or -21). 

There is a class of reference frames which are preferred since the 
laws of mechanics, when expressed in these frames, have a simple form. 
They are called inertial reference frames . In practice no such frame 
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can be found; instead a frame which is sufficiently inertial to 

provide the degree of accuracy required by the problem is selected. 

Examples of sufficiently inertial frames of increasing accuracy are: 

1) a frame fixed on the surface of the earth; 2) a frame fixed at 

the earth’s center with an axis pointed at the sun; 3) a frame fixed 

at the sun’s center with an axis pointed at a nearby star; 4) a 

frame fixed in the ’’fixed stars.” In what follows, a reference to an 

inertial frame will mean a frame selected by a user as sufficiently 

inertial for the problem being solved. Unless otherwise noted, the 

unprimed reference frame will be taken as inertial. 

A physical interpretation of the transformation of velocity between 

reference frames may be obtained by considering just the components in 

physical space. Time is then viewed as a parameter and is not explicitly 

in the manifold structure. Physical space is a submanifold of event 

space, and the . tangent space to physical space is a submanifold of the 

tangent space to event space. Specifically, if M(t) is the submanifold 

of the material continuum at time t (i.e., the time slice at t) , and 
i 3 

V = a — r , i 15 1,. . .,4 is a tangent to M at m f M(t) , then the 
9x 

associated tangent to M(t) is V = a a (t) - — (t), a = 1,2,3 and will 

3x a 

be called a spatial vector , or the spatial part of V. The spatial 
vectors are generally easier to work with and are more familiar, 
however a contravariant vector (e.g. velocity) does not transform 
properly when treated as spatial, and an artificial means must be 
developed for treating them. 
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/ 


The spatial position relative to the primed origin may be defined 


as 


P' = ^'“(x 1 , X 2 , X 3 ) o n ’" 1 — i — = x'“ — , 


9x 


,a 


9x 


,a 


(1.5-27) 


noting the coordinate expression for 31 ’ (1.5-22). Similarly 


P - Ax 1 , X 2 , J 3 )o If 1 -i- - x a A 
•“ a a .a 

dx 9x 


(1.5-28) 


is the spatial position relative to the unprimed origin. Adopting the 
convention of using the same symbol for spatially related fields, 


P' = A*(P’) = p'^x^x 2 , x 3 ) o IT'" 1 o A" 1 — - q - A) ' o A" 1 9 


■ a 


9x 


8 


1 -2 —3. T.,-1 A -1 ,3 A -1 9 

= p (x , x , x ) o IP oA C ,H O.A a 

a 3x e 


f ot .-i ,e A -i 9 
= x o A c o A — 

9x 


a 


8 ‘ 


(1.5-29) 


Recalling the coordinate expression for A (1.5-11) yields 

rj , _ ,..e . A ~ls 9 „ _8 A -1 9 

F-(x - r oA ) — q - P - r o A — ~r 

3x 6 " Sx 13 


= P - R , 

where R is defined by 


(1.5-30) 


R - r 6 o A 1 — . (1.5-31) 

The spatial field R may be recognized as the position of the primed 
origin relative to the unprimed. In summary 


P = R + P’ . 


(1.5-32) 


As noted, velocity does not transform properly as a spatial vector 
field. To circumvent this, differentiation with respect to a reference 
_frame will be defined. Differentiation is conceived of as acting on the 
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cx 3 3 9 ct 9 

vector field as a whole; so if U = u — — • , then — (U) = — (u ) — ~r 

— a dt 9t ^ cl 


Bx 


9t 


9x v 


ck. 9 9 

+ u — ( ) ; The derivative of a basis in its own frame is 

3t 3*“ 

defined to be zero (e.g. |— (-—) = 0 - — - ( — —■) but |— {-?— ~ ) I* 0) . 

3t 3x a 3t 3x ,a at 3x ,a 

With this notation, the definitions of velocity become 

9 


, . a 1 —2 ~3v tt-1 

Vs L (p) = 9P i * > * ) o n 
— 9 1 9 1 


9x 


a 


(1. 5-33) 


and 


v . - 11 ( p.V = V p' V. » 2 ,.? 3 .)„o-I!1 ~. 1 

- 9t K " } 3t 

-1 


ax' 


a 


Inserting the A - related basis, 

-i / a \ 


A 


* U“ 


3 <* ,Y ° Lh o A — — 
9x a Bx 


y 9 

,Y a 9x r 


(1.5-34) 


(1.5-35) 


into the spatial part of the velocity (1.5-26), yields 


9r 


a 


V - 

~ 9t 9x a 


9 , y ,a ,3 

— + c ' c ' v ’ 


a 3 


3x f 


9c ,a 

9 . Y ' 3 *3 

+ c r - X 

Y a 3t 


9x f 


at 


<>y +' r + <£ ^ (c“) *’ 


9x' 


Y 


(1.5-36) 


where it has been noted that 


' a J = 6 y 
3 a °3 


c 1 „ c r = 61 and V T = A. (V 1 ) 


The term c Y t— c’ ■ c Y t*— c p is the component matrix of a "skew- 
ex dt 3 a dt a 

symmetric" tensor (A. 7) since 


C Y <L_ ( 3^ iL 

a 9t v V Bt 


(J J) - f- LA c 6 

\ a a I dt y cl I a 

(6 Y ^) - c^ — (c Y ) 

Bt v } cl dt \ clJ 

- c 6 -2- (cA . 

a at a J 


(1.5-37) 
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A skew-symmetric tensor of type (2,0) is naturally isomorphic to one 
of type (1,0) (in a three dimensional vector space). If this isomor- 
phism is denoted by the "Hodge star operator" (A. 9), then the asso- 
ciated contravariant vector field is 


J i_ ( c 6\ _J 
•a 3t V a/ 3 , 




3x 


,3 


- aj 


,6 


3x 


,6 


* - Q' . 


(1.5-38) 


This vector field may be phy'sically interpreted as, and will be called, 
the angular velocity of the primed frame relative to the unprimed 


frame. The components, Lu 1 , are determined by 

15 y 8 
2 ya c 0 3> 


,6 16 y 3 / a\ 

“ = - - E — \c & ) , 


(1.5-39) 


where E is the "Levi-Civita epsilon” (A. 9). Conversely, a tensor of 

ya 

type (1,0) is isomorphic to a skew symmetric tensor of type (2,0), 
customarily denoted by a tilde, over the symbol for the vector 

field. Thus 


* 





> 







(1.5-40) 


where 


or 


a3 ,6 

- E r CO* 

0 


The double underline denotes a tensor of type (2,0). Combining the 
two definitions yields 

E 


= - e 01 * 5 

Y 


1 y 

- *=- e \ c 

2 e6 n 


fc «)' ' 


(1.5-41) 


Noting the identity, 


a£ y 
e e * 
y e6 


■ ( 5 “ - 6 “ ’ 


(1.5-42) 
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(1.5-41) becomes 

<-,a3 1 f 

“ = y [ 


& a (5® - fi“ 6 B 
e o 6 e 


1 ['■) h ( c 5)l 


(1.5-43) 


Recalling the skew- symmetry , (1.5-37), 


%a3 

or 


cx 3 


«) - - 


£<* 03’ Y 

Y 


TV 3t 

in accord with the definition (1.5-40). 

With this notation, the previously obtained expression for 

velocity, (1.5-36) , becomes 

V = f- <R) + V' + u'“ e x' 6 — — 

■ 3t • “ . ax-“ 

.a3 , . tY 3 


(1.5-44) 


= KT (R) + V’ - £ M UT ' X 1 
3t — — Y 


3x 


,a 


(1.5-45) 


The last term may be recognized as the "cross product" (A. 9) of the 
vectors £2' and Pf at any point in the field. So 

3 


v = ~ r + v 1 + n* x p 1 1 . 


(1.5-46) 


This may be written in the more suggestive form 

ff (P) - |f (R) + f~ (P’) + ft’ x P* . (1.5-47) 

On the other hand, (1.5-32) implies 

ft ( -> ' ft ( -> + ft ' (1.5-48) 

Thus the correct transform for velocity is obtained by setting 

ff <!') - If (£')+£' x P* . (1.5-49) 


Actually no redefinition is required since (1.5-49) holds for any 

i a 3 


vector field. To demonstrate this, let U f = u’ 


of the A-related basis, 


Then in terms 


3x' 
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St - 


at 


- ~ (u ,a ) 
at 


3 

+ u 

,a 

S 

/ 3 


,a 


St 

\Sx ,a 

; 

3 

+ u 

T 01 

3 

/ c » 3 , 

s 

t a 


St 

V 

Sx 1 

t a 

S_ 



3 



St 

( c 

a) 

3x 6 



Returning to the primed basis yields 

fr <u') - <u’) + u’“ c y R | - (c>®) — 

at — at — 3 at \ a/ . , 

ax 


= |1 (u .) + Z't a u > a _ L- 

3t - 3 x' Y 

= (U 1 ) + S3’ x U’ , 

recalling the previous definitions. 

In summary, the proper transformation of velocity, 


V = 7rr(R) + V’ + ft' x P' 
at - — — — 


(1.5-50) 


(1.5-51) 


(1.5-52) 


or 


h ® - if® + if ®> + a' x i ’ . 

is obtained by the definition of differentiation, 


(1.5-53) 


fj- (£') - <U') + «’ X U‘ . 


(1.5-54) 


1.6 Newton's Laws 

Newton’s laws form the basis for much of structural mechanics. 
Although originally stated for a particle, the first and second laws 
have been generalized to a continuum. For any point, m, in a material 
continuum, M, which is in an inertial reference frame, E, 
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( 1 . 6 - 1 ) 


p a - f = a , 

where F is a spatially contravariant body force density representing 
the sum of all long and short range forces at the point; p is the mass 
density; and A is the acceleration . 

The acceleration will be defined as the spatially contravariant 
vector 


A = i. (v) 1L_ 

- 3t . 2 


a 1 — 2 — 3 v 1 TT -1 9 

p (x , x , x )] o H — - 

9t“ ‘ 9x a 


(1.6-2) 


The acceleration in the primed reference frame, 

.2 


A* = ~ (V* ) « - — 
- ^ 3t 2 


I Oty 1 2 3 X ~ 

p (x , X , X ) o II 


.-1 l 

9x' 


a 


(1.6-3) 


is related to A by 


A 1 ^ 

—5 

- 9t 


(R) + V* + Q* x P' 
9t — — — — 


- — 2 <£) + fj <■!') + k' x v' + — (a’) x f' + n' x ■— (p') 


+ 1 •' x (n* x p') 


9 2 ® + + 2 - x I* + f“(^’) x P* + ft’ x (ft* x P') 


9t 


(1.6-4) 


Since Newton's laws hold at any event in a material continuum, 
talcing the "inner product" (A. 7) with arbitrary spatial vector field, 
5P , and "integrating" (A. 10) does not affect the identity, 


fr fv (t) * <( * • -ftfkt) • il) = 0 • (1-6-5) 


The integration may be considered as an integral over a time slice 
followed by integration with respect to the time parameter, or by 
returning to a four dimensional model, as the integral over the 
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material continuum expressed in iterated form in accord with 
Fubini's theorem. This distinction will not be pursued since the 
"associated Riemann integrals", (A. 10), are essentially the same, 

Jv T wfv s Wf» (pa 6p ) 0 U dw 3 dt 

s (M(t)) (fQ 6pa> 0 W " 1 iv * dt = 0 • a ' 6 ' 6) 

The first term is transformed to material event space to remove 
the time dependence of the spatial region of integration. 


&p(%) ° 11 0 d_1 det 


y g (M) 


y _1 detl^— 
U..Y 


dp 3 dt . 


(1.6-7) 


S v 

The Jacobian determinant, det Op /9u T ), may be interpreted as the 
ratio of an Infinitesimal volume in M to one in M. Since the mass 
density, p, may be interpreted as the ratio of the infinitesimal mass 

to its volume in M, the product poll det(^—\ may be interpreted as the 

\3u Y / 

mass density in M and will be denoted "p" . The first term in (1.6-5) 


becomes 


hh 


(p A * 6P) . 


( 1 . 6 - 8 ) 


Recalling that. A * — V , noting the identity 


f^-(V> • fiP - §£ (V • 5P) - v • (6P) (1.6-9) 

and performing the temporal integration (since the spatial region is 
now independent of time) yields for the first term of (1.6-5) 

Ut/s * (S v • <5P) -/ f /g *[$V • fj- (6P)j, (1.6-10) 

It 

where 3T may be taken as meaning . Recalling V = 3/9t(P) , and 

! -i 
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noting that 

2 [h © • h <©] ■ If <£ + &>■ If ^ + i£> 

■ It (-)' ft © + 0( l li£i P> (1.6-11) 

the last term of (1.6-10) may be recognized as the "variation of a 
functional" (A. 11). The first term of (1.6-5) may be written as 

l 3 l/s * ($ V • 8P) - * (| 2 v • v) 

+ 6 /t/s * [t • (1.6-12) 

The term in the inner parentheses of the last integral is constrained 
to be zero, so the term !L may be thought of as Lagrange multipliers for 

the constraint. The term P V_ * V_ is called the kinetic energy density 

A 

and will be denoted 'T". The process which yields the second two terms 
of (1.6-12) from the last term of (1.6-10) may be recognized as a 
Legendre transformation. 

The stationary conditions for the terms in (1.6-12) include 
requiring the coefficients of the independent variables, 6V. and 6_L, be 
equal to zero (e.g. V_ = “ _P and L^ - pv) . Eliminating L as an inde- 
pendent variable by substituting yields 

Us * (P V • 6P) - * (I £ V • v) 

+ 6 A/s + [?- • _(i-|t£)] a-6-13) 

for the first term in (1.6-5). 

The previous results are associated with Newton's first and second 
laws. Newton's third law states that the force on matter, due to its 
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interaction with some other matter, is equal but opposite in sign to 
the force on the other matter due to the first. If the force on the 
matter at tn^ due to interaction with matter at is denoted , 

then 

^(m 2 ,m 1 ) - - £( 11 ^, m^) . (1.6-14) 


At any point, m, in a spatial region, R, let the sum of all long 
range forces be denoted "X", let the sum of all short range forces on 
the matter at the point due to matter outside the region be denoted 
"S M , then the total force, ]?, at the point has components 

f a - x a + s a + y^( R ) 4> a (m, P 1 ) dxdydz . (1.6-15) 


Integrating these component fields with respect to the usual coordi- 
nates, find 



.a 


f op dxdydz 


■/ 


a 


P(R) 


p "^dxdydz 


J\ J(R) 


k a 


p "^dxdydz 


+ /( R )^ (R )^(P 1 »P 1 ) dxdydz dxdydz . 

(1.6-16) 


The last integral vanishes; and the integral of short range external 
forces is customarily replaced by the associated boundary integral, 
yielding 

fu(. R) f0t 0 U " ldxdydZ ”/ M (R) x01 0 y_1 dxdydz + /y(3R) sa dul(3u2 > 

(1.6-17) 

or on the manifold 


* ( £ ■ i) ■/« * i) */-* ( i ■ i ) ■ a - 6 - i8 > 
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1.7 Stress 


In the previous section (1.6) the integrals of the components of 
the total force and the body force vector fields over a region were 
related to the integral of the components of the surface force over the 
boundary. In an appendix (B.l) this result is applied to a tetrahedron 
with an apex at an arbitrary point, m, base normal to an arbitrary 
vector, eN_, and sides defined by the usual coordinate planes. 

The form of the result suggests the definition of a stress tensor 
field of type (2,0), 

v _ 9' „ 3 

= ■ a 75 ® re • u- 7 -!) 


Then a surface force on a surface with an outward pointing normal, 

ct 9 

n — — , is determined as 

9x a 

S = 1 'I , (1.7-2) 

where the "dot product" (A. 7) in the usual coordinates is 


ol a$ 
n a 



(1.7-3) 


It may be noted that the surface force on the face with an outward 

^ 8 3 

pointing normal, — ~ , is - 0 g in accord with the definition 

9x 3x p 

in the appendix. 

It has been shown that there is a tensor field called stress, 
that when dot multiplied by the unit outward pointing normal to a 
surface, yields the surface force density; and that the force density 
so obtained is consistent with (1.6-18). However, this has only been 
shown for the usual coordinates. In Appendix B,2 it is demonstrated 
that the relations are true for arbitrary coordinate systems. Thus 
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the relation, (1.6-18), may be written 


/» *( £ ■ £) ■/* * ( s ■ io 'f* * ( a • 1 • i) ■ a -’-" 

The last integral can be converted to a volume integral by means 
of "Stokes' theorem" (A. 10), 

/ 

3 




d * {N • JL 


•r. a 
dx 


*d* ( N * £ 


3x 


= /■ * 

,/R 




a 


(1.7-5) 


Since the region of integration is arbitrary, the components, and thus 
the vector fields, are related by 

I - X + V • £ , (1.7-6) 
where "nabla, V," with a dot denotes the "divergence" (A. 12) . In the 
usual coordinates, this is 


M - x a J- + JL (a 3a ) JL 

r dx a dx a Bx^ 3x a 


(1.7-7) 


1.8 Virtual Work 

Substituting (1.6-13) for the first term and (1.7-6) for the 
second term of (1.6-5), yields as a consequence of Newton's Laws 

- * & + S fifs * [pI • (l - h *)] 

+ I 3 t/s * 

-Jt/hM * * [(v • l)- fip) = 0 . (1.8-1) 

Body forces generally depend on the amount of matter present, so 
it is convenient to transform to material event space. The product 
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X 


Bet IZ£\ 


\ g 3 / ±S called the body fQrce referred to the undeformed 


volume and is denoted "X. " The body force integral may then be written 

-/t/hm * & - s ?> - -fjfs * <2 • ■ a. 8-2) 


The body force may be separated into two parts, a conservative part 
which maybe obtained as a "gradient" (A. 12), C x = - V W , and a non- 
conservative part, N X. The body force integral (1.8-2) becomes 

~Jt fu(t) * ( - • “ -fx fi * (- £W' • 6P) * ( N X • <5P) . 

(1.8-3) 

Noting that 

- • 6P - - W(F + 6P) + W(p) + o(| |«p||), (1.8-4) 

the body force integral may be written 

t) * - ' 6 - = 6 fijs * 'frfs * (N £ ' SD- U-8-5) 

A 

The function W is called the potential energy density field. This 
result may be generalized slightly by allowing S to be a non- 
conservative explicit function of time, in which case W is said to be 
lamellar or irrotational . 

The integral involving the stress may be rewritten by noting the 
identity 


(i • JL ) • - M = V I • sA - s : V 6P . (1.8-6) 

The double dot means contractions are to be performed on the first 
indicies of both tensors and then on the second. In terms of the usual 
coordinates this identity has the form 


a 3 


3 


- -yr Co®* «P 3 ) - ® 


3x 


a 


ag 3 


3x 


3x 


a 


(Sp 6 ) 


(1.8-7) 
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The stress integral in (1.8-1) may thus be written in two parts, 

-/f/ M (t) *\(z • £) • & z] ■ 't 1 • 6 ^)] 

= V 6p) . (1.8-8) 

The first stress integral in (1.8-8) may be converted to a boundary 
integral via the generalized Stokes' theorem (A. 10), yielding 


-/f/ M (t) «*)] = -frfm(t) * (**£•«) 


. (1.8-9) 


The boundary integral may be separated into two parts, the spatial 

p 

surface, 3M(t)g, over which surface forces are prescribed, and the 

p 

surface, 3M(t) , over which positions are prescribed, P. The condition 

r 

that actual positions equal the prescribed can be introduced with S_ as 
a Lagrange multiplier. With these definitions (1.8-9) may be written 

-fr/ni t)* [l ■ (£ • »)] “ • 6 1) 

(s. • • (t - p ^)] • (1 - 8 - 10 ^ 

The second stress integral in (1.8-8) may be transformed to 
material event space yielding 

(l : i -fijj *[(^ : z & *)° n l det (^) l] •< 1 - 8 - 11 > 


The stress tensor may be expressed in terms of the usual basis for 
material event space. However, the stress represents the magnitude of 
the force on a unit area perpendicular to a unit normal in S, not M(t) . 
The rectangular parallelepiped with a base defining a unit area and 
with a height determined by the unit normal becomes, in S, an oblique 
parallelepiped with volume equal to the Jacobian determinant. To 
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obtain the force for a unit normal and unit area in S, ^ can be 

multiplied by the Jacobian determinant ! det i* The result is 

called the stress referred to the dimensions of an element before 

A 

deformation, and is denoted " E. n . 


The differentiation in the gradient may be converted to material 
event space coordinates. In terms of the usual coordinates, this 
becomes 
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. (1.8-12) 


Noting the identity 

; 86 


r, (x^olloy \>yoII ^oy 

6 „ ± \ K . = _°L 


3x 3 du* 


OjJ 
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(y 6 oil ^oy ■*■) a (x^olloy *) __ 

TT ^ 


3u 
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oy 


oyoll 


3u 


a (?oH X ) (x 6 oIl) 

6 -2- oh -1 


3x^ 


3x 


(1.8-13) 


(1.8-12) may be written in the form 
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O n o“ e o n iii 


9 (x o IT 


=H)- 


(1.8-14) 


where the first term in the product may be recognized as the related 

— AV£ 

stress in S, a 1 . Alternatively, the expression for the related stress 


in M(t) , 


n* p P ® ~) - 

\ ax' 8x / 


- o ir 1 ] n -ifi^nio rT 1 ]-^- ® , u.s-is) 

ax' J 9 x £ J 3 x 3 x p 

could be directly substituted in (1.8-12) with the same result. 

The second term in the product in (1.8-14) may be recognized as 
the variation of the strains . The strains are determined from the 
fundamental bilinear form (1.5-3). The Jl-related form in s' is 


$ s<* a ° i( * S ° n.) d ^r d ^5 

aS s-y ax w dx • 

dX* dx 


(1.8-16) 


The strain tensor , n E", is defined to be 1/2 the difference between the 


E - i [n*(b) - b] 


i |a(x a o ip 3(x a o ip 


- 6 d^ ® dx 5 


(1.8-17) 


the usual form can be obtained by noting 


5 Ot r,—CL 

? = ox 

Y -—6 

3 x ox 


(1.8-18) 
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The strain can then be written in the form 




3(x a oH - ?) 3 ( 

a „ — ol 

X oh - x ) 

■ 

3^ 

,6 3(x a 0 n - 7*) 
a 3^ 

d^® 

[3(x Y oII - ) , 3(x^o]T - x^) 

H 

bT' 

3(x a olI - x*) 3(x°oIt - x*)] 


a 




a. 


-6 

3x 


rr 


3x 


a? 5 




(1.8-19) 


-6 


In the customary notation, the difference in coordinates, x oil - x , 


6 X 

is called the displacement , u , the form — 


3u | 3u Y 
3 x^ 


is just the 


usual definition of the linear strains, and the term ~ 


~ a ^ a 
3u du 

37 s 3^ 


is 


just the correction to provide the geometrically nonlinear definition 

< 

of the strain. Now* noting that 


3(p^o u) 3(6p^o y ) + 3(5p^o y) 3(p^o y) 


3^ 


3^ 


3x 
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3(p 6 + 6p 6 o .y) 3,( p 6 + 6p 6 o y) 3(p 6 o y) Mp 6 Q W) , C (| |^ r | |) 

3? 


3? 

6. 


3^ 




6 


J 3(p o y) 3 (p o y) _ & 

ey 
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6 2e 
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( 1 . 8 - 20 ) 


(1.8-14) may be written as 


f . . 


AY£ + Jey| I 3(P 6 Q U) 


3x 


3(5p^o u) 


rr 


3x 


( 1 . 8 - 21 ) 
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Since the stress tensor is symmetric, interchanging indices gives 



3(p y) 

8(6p^o y) 9(p 5 o y) 3(5 p^o y) 

15 J 

[ 3x 

■ 3 ? ; ' 1 3^ fi? J 


= i$ Y£ 6[2e 

2 ey 


] = S Y£ 6b 


ye 


( 1 . 8 - 22 ) 


Thus the second stress integral (1.8-11) may be written 


fifn( t)*(l : I & l) =fcfs : 6 i) + - l)]> 


(1.8-23) 


where the term in the inner parentheses of the last term is constrained 

A 

to be zero, so the term £ may be thought of as a Lagrange multiplier. 
The bar over nabla implies the associated expression in material event 
space. 

The stresses may be represented as the sum of non-conservative and 
conservative parts. If the deformation is adiabatic or isothermal, 

then the process is conservative and a Y£ de is a perfect differ- 

y£ 

ential (Washizu, 1968). This assures the existence of a strain energy 
density function A. The second stress integral (1.8-23) may thus be 
written as 

f f *( Z : V6p\ -&{_(_ *(b + 6 f_ f *rl ! (V P - E)1 

•'T •'M(t) \ ' • / T* / S * / T* / S L “J 

+ f?fs * ( N ” : 5 S) > (1.8-24) 

where the prescript "N" denotes the non-conservative part of the stress. 

In summary, Newton's laws have as a consequence (1.8-1), which 
may be rewritten, inserting Equations (1. 8-5 , -8 ,-10 and -24), as 
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- 5 /?/s * & + ! hh * <«> + hh * ^ 

+ «/ Y / ¥ *\h • (l-h4 + 6 /t/s *[^ : (I I - I)] 

--^/M( t ) p *j 6 ^ *(t - P 4]| + I m/? *( fjV - ' 6 *) 

-/t/shm *k ■ 6 

b p 

-fxfs *( N * • 5 ?) + hh * ( N * : 6 i) -■ 0 • 


In the absence of inertial terms and non-linearities, and assuming the 
prescribed surface forces and displacements are constant, this becomes 

*f= 


fi j/s *"<*> + fs * <«>. + /s * [i : (^ t - I)] 

h {jfe )p * •(*. - P 41 - jL( t ) s * 6 ( p ^ • *) 


= 0, (1.8-26) 


which may be recognized as the Hu-Washizu variational form. 

1.9 Virtual Work for an Arbitrary Reference Frame 

The expression of the principle of virtual work in Equation 
(1,8-25) was based on Newton's laws in an inertial reference frame 
(1.6-1). The result may be transformed to an arbitrary reference frame 
via the results of Section 1.5 (Equations (1.5-32, -46 and -54)). It 
may be noted that since the strain is defined in terms of bilinear 
forms, and since bilinear forms are independent of reference frames, 
the strains do not change. Proceeding term by term. 
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■ /t /? *[l R + I' + Q' X P') .(Ij r + V- + a- 



■ ft* • + 1’- 1'+ ia'x p'l-ta'x p’]) 

+ p{ft^’ 1' + I!f *] • t«’ *£'1+[£'H8' X P'j) . (1.9-1) 

Noting the identities in Appendix A. 13, this may be written 

frfs H ^ = -/i|i at - I i-Js * ( ^ + l/s ' D 

+ Js *(p [(l ,# I')s “ l r L’]) • E f 

+ h - fs * ( ^ ?) + h ± x ^'fs * ( ^’> + ^ -fs *(PZ Tx V’)} , 


(1.9-2) 


where U is an identity tensor (i.e., u • x ■ x) . The potential energy 
becomes 

Jjfs *[W(P)] =/y/s * [W(R + P’)J 

=Jt/s * £.')]• (1.9-3) 

The strain energy, as previously noted is independent of frame. The 
constraint between velocity and position becomes 
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fifvfc' (*‘M] 

“A/H" ft - + r + “' **')• (fea + s 

~|t - a' xr)J 


f + ft* x P' 



Tl/'l*^ 


'• i') + %£-fs * ( p- ,) + s’ • *(pp' >< v') 


- **• l?7? * < ^ ,) -/? *(&•• |j a-) - 


The next Integral becomes 
A 


hh * [ ' < vp - S)1 -frfs * : - l) • 

J 


(1.9-4) 


(1.9-5) 


since _R is constant spatially, and since the stress and strain are 
independent of reference frame. 

The integral on the temporal boundary becomes 


yf/s *(PV- 


<5P) = 


r f§ *[p(ft R + v- + a- x p-) • (sr + Sp'j] 

“ larjft x «•> 

+ |f R • Jj *(pfip’) - a' X SR- yi *(pP') 

+ 6R* /”g *(p T)+Jj *(p V 1 • SP’) . (1.9-6) 

The prescribed positions on the spatial boundary become 

-hfmtAi-k - P l)] ■ *[£'•(« *f-l- V)| 


-fU 


T(/3M(t). 


* (S' ‘ P f ) 



The other remains unchanged. The prescribed forces on the spatial 
boundary become 



(1.9-8) 


The non-conservative body forces become 

-hh ~hh * + «t’>] 

- - Js *(”!') +_/§■ *( N 1' - «£') • (1 . 9-9) 

The non-conservative stress integral remains unchanged. 

The dot products of R with vectors defined in the primed basis 
introduces the direction cosines into the variational form. The 
relationship between the direction cosines and the angular velocity 
(1.5-39) may be introduced as a constraint 

5 /f * 

The Lagrange multiplier, H ’ , may be interpreted as the angular 
momentum. 

Inserting Equations (1.9-2) through Equations (1.9-10), yields 
for Equation (1.8-25) 


[h-. E 5 c Yl_c“) 3 1 

p 2 ya 8 at 3/^, fir 


(1.9-10) 
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SECTION 2 


SUBSTRUCTURES 


A substructure will be defined as a finite collection of material 
continua satisfying certain properties. The continua must all coexist 
in the same region of time. They must have no interior points in com- 
mon, however, each must have part of its boundary in common with 
another. The shared boundaries of the continua will be called the 
internal boundaries of the substructure, and any parts of the boundaries 
of the continua which are not shared will be called external boundaries 
of the substructure. The union of the closures of all the continua 
must be connected, and the boundary of the union is the union of all 
of the external boundaries. The position and stress are generally 
required to be continuous at all internal boundaries. 

The final requirement is that the positions deviate from a 
reference spatial state by a small amount for some reference frame. 

This last requirement forms the basis for an approximate solution, 
which is the subject of this section. 

2.1 Finite Elements 

The variational expression (1.9-11) has as known functions of 

the coordinates and the other variables the potential energy, W f , the 

. A , NA 

strain energy. A, the non-conservative body forces and stresses, X 

N A 

and JL , and on the external boundaries, the prescribed positions, 

P P 

and surface forces, P^ and S_. The prescribed positions and surface 

P P 

forces on the internal boundaries, P and S_, the velocity, V 1 , the 

a ; 

position, P ’ , the stress, JL , the strain, E , and the surface force, 
are all unknown functions 'of the coordinates. 


““•*^*91 *ut 
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Direct solution of the governing equations is seldom possible. 
Instead, the unknowns are represented as a series, and the least sig- 
nificant terms are truncated. The domain of the expansion may be the 
whole substructure, the material continuum or a mixed form involving 
both. Moreover, the domains for the various fields need not be the 
same. When the domain of the expansion is the material continuum, it 
is called a finite element . Since use of expansions for the other 
domains is rare, only finite elements will be treated. 

It may be noted that the series expansion could be in the 
temporal as well as spatial coordinates. In fact this has the advan- 
tage of reducing the partial differential equations to algebraic 
equations, as opposed to the ordinary differential equations 
obtained by expanding only in the spatial coordinates. However, the 
temporal approach will not be pursued here, and the conventional 
spatial expansion will be used. 

The development of a finite element seems, at present, to be as 
much art as science. The developer may arbitrarily select the geometry 
of the elements, the fields to be approximated, the approximations 
to be used, and the variational form to be applied. Aside from certain 
mathematical and practical considerations, the worth of an element is 
determined by the accuracy it provides for a given number of degrees of 
freedom, and by its convenience in use. Although the development of 
finite elements lies outside the scope of this presentation, some of 
the more pertinent aspects of the method will be included. 


46 



Geometrically regular material continua are generally selected 
because of their analytical convenience, and the ease with which they 
are assembled. Parametric elements are a notable exception, where an 
irregular spatial region is mapped into a regular region in a parametric 
space. 

The choice of which fields to approximate can be divided into 
three categories: first assumed compatible displacement fields, 

second assumed equilibrium stress fields, and third mixed partial 
assumptions on the stresses and displacements. There are further 
choices within each category. For instance, with a displacement 
assumed model, the strain could either be separately assumed, or 
determined from the strain displacement relations. In fact, the dis- 
placement fields could even be different depending on application 
(e.g., one for strains, another for the potential energy associated 
with body forces) . 

The fact that any complete series expansion provides an exact 
solution to the equations, makes them desirable choices. (However, 
any series expansion may be acceptable.) Polynomials are frequently 
used because of their convenience and since (the set of all poly- 
nomials to order °°) is a complete function space for the class of 
continuous functions (Weierstrass theorem). Trigonometric expansions 
are also used, especially for the circumferential direction in 
axisymmetric elements . 

The variational form used may be some variant of (1.9-11), 
which is essentially of the Hu-Washizu form, or some variant of the 
dual form obtained by a Legendre transformation of the strain energy, 
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yielding a complementary energy formulation. Generally, finite 
elements are developed ignoring the kinetic energy terms, considering 
only the strain energy and sometimes the potential energy. When 
the finite element developers do include kinetic energy, it is 
usually in the simpler form of Equation (1.8-25) rather than Equation 
(1.9-11), with the velocity field explicitly assumed to be the time 
derivative of the position field. In this special case, the associated 
terms are usually referred to as "consistent mass." 

To be of practical value, the finite element field approximations 
must converge to the correct values for the field as the element size 
decreases. A finite element generally must satisfy the following 
criteria to be convergent. First, the approximation should be capable 
of representing a constant value of the field or its derivatives up to 
the highest order appearing in the functional. Secondly, the approxi- 
mation should be continuous to one order less than the highest deriva- 
tive appearing in the function. Generally this is met by requiring 
that the highest order derivative be continuous almost everywhere 
(i.e., except on a set of measure zero, e,.g., piecewise continuous). 

The finite elements are assembled by applying the variational 
form to the whole substructure. The conditions at the interfaces are 
just the previously mentioned continuity requirements. If these are 
met, the internal boundary integrals cancel. However, the requirement 
may be relaxed, If the discontinuity is at the boundary, by explicitly 
retaining the internal boundary integrals. In practice, when using an 
assumed displacement field, the satisfaction of compatibility require- 
ments is usually made automatic by transforming from coefficients of 
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an approximate expansion to the displacements at nodes on the 
boundary. If the displacement on a boundary surface depends only on 
the displacements of the nodes defining the surface, compatibility 
is established by requiring the displacements of those nodes be 
identical for adjacent elements. This is the basis of the "direct 
stiffness" method. 

As an illustrative example, the derivation of a beam finite 

I 

element based on an assumed displacement field is presented in 
Appendix C.l. The derivation includes the effects of geometric non- 
linearities. 

2.2 Finite Element Dynamics 

In this section, and in the material that follows, it will be 
assumed that the assumed displacement finite element method is 
being used and the assumed displacement field is known. Alternatively, 
a mixed or stress assumed formulation may be used, provided that the 
displacement field is known, or can be obtained in closed form. To 
incorporate dynamics, the variational expression in Equation (1.9-11) 
must be used instead of Equation (1.8-25). Taking the indicated 
variations in Equation (1.9-11) and noting the definitions of vari- 
ations of vectors in Appendix (A. 14), 
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(2.2-1) 


Noting the derivative relationships, 
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A - 6 R = < A - 6 R> 6R , 

A T • <5 T p* = (a’ * 6P’) - A' • (SP 1 (2.2-2) 

and integrating by parts, yields 
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Equation (2.2-3) Continued 
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- ‘j \ 4 It ■=:- lr[ h i t 
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t| lT( h i I -yp 'I *7? *{» It p’“ * P*'V*cJ 

+ p(S' X P’) Y p' S c Y })«c e . 

— aj/ a 


(2 . 2-3) 


The stationary conditions of this form provide, in essence, the Lagrange 
equations for the element. It has been assumed that all vectors except 
R. are expressed in the primed (body fixed) vector basis, and that R is 
expressed in the unprimed (inertially fixed) vector basis. 

If the displacement and velocity fields are assumed to have the 

form 

9 . - 9 -9 


pt _ ,, ^ , 3 , — 9 , ^ i 

- x 3x’ + y + z + £ ±1 j. 1 ! and 

i=l 


N 


- = £ s 1 ’ 1 P 

j=i J J 


(2.2-4) 


some of the stationary conditions may be immediately identified. The 
coefficient of fill must vanish over the time domain, T y 
2 

g * ■ 9 1 ___ 

M r R + M ~ ft 1 x ( P 1 + . P ’ .q) + Mft‘ x P’ q 

— at — * 0“ X— 1 — 1— x * 

+ M ft’ x [ft 7 x ( "P* + .F’.q) ] + M.P* .q + M ft’ x J V 1 q 

— — o— i— i n x— x M — * i— i H 

m fs* |- 1 *' + % + Xt(t) p *^] + ^M(t) *(V) . (2.2-5) 
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where the repeated Latin subscript implies summation over the M q's, 
and 

m =/- * {&}, 0 p- =/- *{a(x + 7 + 7 yj /M and 

J’ =fs * • 

In addition, the coefficient of 6R must vanish on the time boundary, 

3T , 

j“ j P i— (2.2-6) 


where 


jV* * {frjvn/M . 

The coefficients of the independent variations in 6’V’ (i.e., the 
<SjP> must vanish for all time. 


N 


5 ij V i P ki G 


i=l 


( 2 . 2 - 7 ) 


where 


ij 


v = fs * { M' ■ jl' } > ij G = fs . 


The coefficients of <Se 6 2 ,6'_S' , and the terms VW’-6'p', 

N A A p 

21 ^ .L > I : > S/ ' ^'P', ar *d S'-d'P' are the terms associated 

with the elasticity problem. It will now be assumed that the elasticity 
problem has already been solved by the methods noted in Section 2.1. 
Moreover, it will be assumed that the results are of the form 


i Q = i Q + y V + y c j* > 


(2.2-8) 


where the are the generalized forces associated with the surface 
and body forces, and the _.q are the generalized (generally nodal) 
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displacements with any constrained additional coordinates eliminated. 
The coefficients are presumed to include the effects of preload 
(geometric non-linearity), and any other effects, linearized for small 
q, so that they are constant or at most known functions of time. 

j 

If additional displacement degrees of freedom were eliminated 
during the solution of the elastic problem, those constraints are 
modified by dynamics and the static coefficients are improper for the 
dynamics problem. Properly, those extra coordinates should be retained 
in the problem; however, the problem will be disregarded here by 
selecting the field resulting from eliminating the constraints as the 
initially assumed field. In this case, the previously noted terms can 
be replaced 


fi/A 


a“ g + 


a3y<5 A a $ N/\a$| 
a e y6 - o + a ) 


<$£ 

a3 


(VP/- E): 6 f + £ :6V p* 



Returning now to the stationary conditions for Equation (2.2-3), 
the coefficients of the independent variations in 5 1 P/ (i.e., the 6.q) 
must vanish for all time. In the simplest dynamic case, 0/ and R both 
zero, this gives 


i Q = i Q + u K j q + ij c + ik G k p • 


( 2 . 2 - 10 ) 


Solving Equation (2,2-7), assuming _^_.V is non-singular, or taking a 
pseudo -inverse if it is, the stationary conditions are 
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( 2 , 2 - 11 ) 


, .M .q + ,,C ,q + . .K .q = .Q - °Q , 

ij 3 ij j ij 3 i i 


where 


ij 


N N 

m - Z £ 

k=l £=1 


lk G k* V ji G ‘ 


This result is the usual dynamic model. If ,V* = P' then the mass 

J“ 3 “ 

is called consistent. If a piecewise continuous linear field is used 
for the velocity field, the result is nodal masses. However, just the 
assumption of this simplified field does not eliminate coupling between 
the nodes. To eliminate coupling the constraints must be violated. 

This procedure is justified if the resultant energy error vanishes as 
the element size tends to zero. It may be noted that the assumed 
velocity field appears nowhere else. Thus for what follows, the 
question of proper velocity field is immaterial and will not be 
further treated. Instead, it will be assumed that the array M is 
available. 

The mass terms associated with the displacement field could 
conceivably be derived based on a different assumed displacement field 
from that assumed for the elastic solution, with or without violation 
of the constraints. This possibility will not be pursued here. 
Returning to the general case, the stationary conditions of Equation 
(2. 2.-3) provide, instead of Equation (2. 2.-1), as coefficients of 
the independent variations in 6 f P^ (i.e,, the 6.q) 


- ft 1 -. I ,1 • O'.q+M .P'* ~ R + 2ft* • ,0 .q + -r— ft\0 

- - - ji= - 3 i” 9t 2 - “ ji- 3 ~ i- 


3' 


3* 


— ft'*.. 0 .q + . .M .q + C .c[ + K.q = ,Q - .Q , 
3t - jx- ij ij ij i x i v » 


( 2 . 2 - 12 ) 
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il 'fs •■$«(£'- iI')2-o£' i^’ ]} ’ 

ij l =/?* tOUiP'-jP'Ju - £■ jP'3> , 

,0 = f- * {& P’ x .P’} and 
3. j u (j i** TT 

i.j— =/? * * j !' 1 • 


ij 

The remaining terms govern the rotational behavior of the element. 

The coefficients of 6 H* must be zero so 

o 


,6 1 6 y 3 a 

oj = - -? c c! tt- c n . 

2 7 a 3 3t 3 


(2.2-13) 


This may be recognized as the nonlinear differential definition of fi 1 . 

The coefficients of 6'Q' must also vanish, 

H f = ft 1 • I + ft 1 • I q + I • ft’ ,q + Q f • ..I .q ,q 
“ - = - i= i i= — i — xj= i j 

+ M n P , x| 7 R + M.P’ x|~R.q+,0 .q + . .0 4 q d 
0 — dt — jl 9t — i i— 13 - i 


where 


1 K 0 p'. 0 p-) u - 0 p- 0 p-]} . 


(2.2-14) 


The Lagrange multiplier, H' , may be recognized as the angular momentum 
about the primed origin (±,e. 9 J~ *{$ p f x V}). 


o 

The variations, <5 c are not independent; only three are. Instead 
of reducing the generality by selecting a specific three parameter set 
(e.g. a set of Euler angles), the three independent equations can be 

obtained by noting the skew symmetry of . The difference with 

a a 

its transpose will contain only the three independent linear combinations 
of variations. 
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(2.2-15) 


a * 3 a y Y * 3 

a a o c • = a n c 1 cl 6 cl 

3 a 3 a 6 6 


= c Y 6 3 6 C cj 6 c x 

6 a e y 6 6 


a y i ^ £ x 

= - a„ c 6 6 c P o c 


C 

5 ’ 


i-rs v, v_» ^ r 

3 a e y 6 

where the skew- symmetry has been noted. Adding these expressions, 
noting appendix Equation (A. 9-14) and interchanging indices, 
a „ B lT/av nft\ / n c » eV] 

VJ 


B 


5 < ■ i [(*2 4 f) f 


c ~ el 
eC 6 <5 j 


= b 


n 


(t 


n c j? £ 
e £C c 6 5 c 6 


)• 


(2.2-16) 


The quantity in the parentheses represents the three independent 

variations (T) = 1,2,3), so the stationary conditions h 1 "" 1 = 0, can be 

obtained by multiplying the coefficient of 6 c^ by c Y . 

a J y a 

3 

The coefficient of 5 on the time boundary, after multiplication 

. n3 £ , 

by e c is 
e a 


n -Yl3 £ ( u f 16 

0=e c [ h’ ~ 0 ,, c 

£ a \ 6 2 y3 


= (s 6 6^-6 6 

\ e Y £Y 


y + /? 

a JS 

R’ I * fi 


A, 3 d a, Y f Y Y,^, . v Y 

PP r + c T v +c f (^xP)' 

3t a a — — 


+ 6 ey (fi* x p T ) 




r lj 

■f 


V s 


£ 3 a, .£Y ,y 

c aH r + 6 v 


It r “ + V ’ £ + 


;«• x p') e jj 


S’ -fs* !&I' x [IJR + »■ + 0* x P']j , 


(2.2-17) 


which is the same as Equation (2.2-14) (i.e., H 1 is the angular 

momentum about the primed origin on the time boundary as well as the 

interior, t ( T) . The final stationary condition is that the 

coefficient of 6 cf, after multiplication by 0^^ c £ , must vanish in 

ot 0 a 
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the interior of the time domain 


' B c E /-h' i £ 6 J-h> i E 5 8- c Tf 

: o | S 2 yMt a St { 2 yS o 6 2 yB 3t c a 

f *( A 3 a r , ,1 B a 3 a 3 ’ ,8 

J s T p 37 r L- * - j _ P 37 r 3l p 

+ [vfi’]%’ e P - e c 1 } 


P 3 t r 3t p 


+ /M(t) p *t- S ' Y P P' S C a| + / M (t) s *(- Ps ’ Y P’ e <£)) 

— k h ' S j(4 S r - 6 E y + h '6 ( 5 E 5 y - 

+ /s ‘I' I']' Ij- r“ - e f $ |i p.^ ^ r“ 

+ E f P ' e (vft- e - N 5 £ )j 


+ f *L F nB P n* B = 

+ >(t) 0 ( E £ P S 


'SM(t) 


E 716 p' B S ' E ] 
£ 


=(fj- H')^ Cg' x H’) n xj ^ *[- p [(R 1 x P') x Rp- p 
+ [p’ x (VW' - N X’)] n } 

+ V x i'f + fmt) s *{- ? x P «*| • ( 2 - 2 -i8) 

Rearranging and multiplying by -^Tfy gives 

fj- a’ + g’ X H' =J- *jp[j~- P’ + R’ x P'J X R + p’x|^- vft’+ N x'|J 

X S'\x f dH(t) s *\z' x V) . (2.2-19) 

This may be recognized as the identity 

hfl* f p ’ X £v> I^P' x [|- R + p -J + p . X 0a| . 
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where it may be noted that pA = F^. The expression can be converted 
to a function of ft' , instead of JH 1 , by inserting Equation (2.2-14) . 
To linear order in q this is 

4r ft’ • I + ft’ x (I * ft’) + ~- ft* * .1 .q + ft’ • .1 .q + ft’ x (ft’ • I) .q 
dt — = — = — dt — 1= l “ l — — i= 


+ I • ft’ q + .1 * ft' . q + ft’ x .1 • ft’.q + M^P’x - — r- R 

i= dt — i i= — l — — i 1 0— rs 2 ~ 


+ M(ft ’ x C P') x R 


I 

3t‘ 


^ ^ 

+ M ,P’ x ^ R q + M , P ' x R . q + M(ft T x .P’) x ~ R a 

i Q t z — ^ i d t — i — i — d t — i 

+ , 0 q + ft ' x 4 0 . q 
i— i — l—i 

= M P'xf- R .4 + M(ft'x n P’) x R + M(f2'x .P') x R .q 
1— at 1 — U— d t — — i— dt — 1 

+ fs * {£’ X [- vO- + N I-]} + f mt *{ P P'* S'} + f mt) *{p' X V}. 

P s 

(2.2-20) 

or rearranging, 

If R’* I+5'xI • ft’+ |i «’• ± I ^ + £'• ± I i k+ g' x (£'• I) q 

+ iJ ’ ft 


+ I • O' q + O’x I • a’ .q + M P' x R f M ,P’ x R q 
1 - i- - 1 0- g t 2 - l- 3t 2 - i 

+ .0 q + ft’ x .0 . q 
1“ i — l— 

= fi * { £' X t- V&'+ N |-]} +f 3m) * t P P’x S'} + J m(t) *{p'x P S'} . 

P S 

( 2 . 2 - 21 ) 
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2.3 Finite Element Assembly 

A substructure, as a collection of material continua, is itself 
a material continuum. Therefore it must satisfy the same variational 
principle as the finite elements. Each finite element deformation field 
may be thought of as a piecewise continuous deformation field for the 
substructure, defined on the element and zero elsewhere. The varia- 
tional form for the entire substructure differs from the sum of the 
forms for the finite elements by the presence of pairs of internal 
boundary integrals in the later. The displacements must be continuous 
(compatible), and the stresses must be continuous (in equilibrium) 
throughout the substructre. In particular these conditions must be 
satisfied on the element boundaries. In this case, the pairs of 
integrals cancel. Since the finite element method is only approximate, 
there are generally violations (usually of equilibrium for the dis- 
placement method) and the integrals do not vanish. However, under 
appropriate conditions, it can be shown that the contribution of these 
integrals goes to zero as the mesh size approaches zero (i.e., as the 
number of elements in the model approaches infinity). The details of 
these effects belong to the study of finite elements, and will not 
be treated here* It will be assumed that an appropriate model has been 
developed, and that the mesh (number of elements) is sufficiently 
refined so that the contribution of the integrals is small enough 
to be ignored. 

The result of the application of the compatibility requirements 
to the collection of finite elements is to assemble them into the 
substructure by eliminating the dependent generalized displacements. 
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In the case of nodal displacements, the assembly of the stiffness 
terms is quite automatic and is called the direct stiffness method, 
lo generalize this procedure slightly, all the parameters associated 
with a given finite element will be denoted by a pre— superscript index 
associated with the element. For an element, a, the vector ^P T will 
be divided into two parts, a constant portion, ' , defining the 
position of a local origin (used in developing element properties) 
relative to the substructure, plus the position relative to the local 
origin for the element, P 1 . Now if the variational forms (e.g. 
(2*2—3)) are summed for all of the elements, assuming the internal 
boundary integrals are negligible, the stationary conditions for the 
frame variables are the sums of the element conditions. If the 
compatibility requirements between elements are expressed in terms 
of N independent ^q in the form 

N 

ot a a 

i q = ’ < 2 - 3 - A > 

they may be incorporated in the form with Lagrange multipliers, 

E Mr* 

« s « ^ CX. >0C Ot As 

where E is the number of elements. (In the case of direct stiffness, 

the q are the nodal displacements in global coordinates, the ^q are 
3 i 

the nodal displacements of the a element in local coordinates, and a T 

ij 

consists of direction cosines.) The stationary conditions on the 
deformation variables are the same as for the elements individually, 
except for the additional term The coefficients of 5^A just give 
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Equation (2.3-1), and the coefficients of the variations of the global 

coordinates, 6 are just ?.T. 

J a=l i=L 1 lj 

Summing Equation (2*2-5) over all elements provides the trans- 
lational' equation, 




3t' 


— + fr 2.' x o-' + * & £' + 2 «' * •!' j4 

+ n*x(o' x s Q l') + a; x (a-x ®p') ^ + ®p> = s £ , (2.3-3) 

where double Latin subscripts imply summation i=l,.*,N and 

E 

s. 


5 M = £ a M , 
a=l 
E ’ 


0“ \S0- E0- I 


'"’V S M , 


E E “M “jI/ S M and 
a=li=l 


37 M 

E a 


s i - E(X- * { - ^’ + ^ + f 1 A *{VA , 

a= ‘ 1 '* / s ^ M(t)p h\(t) 5 ! 

noting, however, that only the external boundary integrals are 
retained in the force calculation. Similarly, the sum of the nominal 
rotational Equations (2.2—21) over all elements gives 

S i • fj «’+ S’ x s i • n- + II n- s i A + ?i • |1 a- A 

oz - - at - j= y at - j H 

+ £’ • ji + jj *«’ jQ + £' x(g' • jl) ^ + Q' X *1 • a* 4 

2 2 

+ S M ?P^ i— p + S M %*x JL- r X + ^0 .§ + ft* X S 0 4 - S T 

J 9t J J J “ J" J “ 

(2.3-4) 


M or x ^ r + 
3t 


where 
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again noting that only external boundary integrals are retained in the 
torque calculation. 


As previously noted, the incorporation of the constraints via 

Lagrange multipliers results in the appearance of the additional term 
ol . th th 

j_A in the i equation for the a element. Solving these equations 



E M 

a 

L r (Y Q 

2-r jA T = 0. In other words, the new stationary conditions are 
a=l 1=1 1 

ot 

obtained by premultiplying by T and summing over all i in all 
elements a . The elastic equations (2.2-12) become 
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( 2 . 3 - 5 ) 


£'• j I • O' - «*• 1 3 1 • a* , 
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3t 
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+ 2 


» . 
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- *j£ ^ + ^ + ^2 • 
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3t 


9 ot ^ 

£ »q 


+ iSL 




where 






S 

j 


Q 



again including only external boundaries. 

It may be seen from the above equations that the following 
parameters are required for each finite element: 
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a p' “pt a T a T a 0 a T a n Sr a r °v aa n 
* EO - » i- ’ i ’ ’ i- > ik= ’ ik£ * i!? 1 ’ ik C ’ ik*’ i Q ‘ 

1 (Y OtO 

The last two parameters, ik K and Q t were calculated for a beam 

element in Appendix C.l. The remaining terms are calculated for the 

same element in Appendix C.2, as an illustrative example. The 
ct ct 

parameters, ggP f , appearing in the preceding equation are 

specifically associated with the assembly procedure. It should be 
noted that the positions are defined for a deformed state which may 
require a geometrically nonlinear static analysis for their computation. 
Finally, the generalized force parameters must be computed from an 
assumed environmental model. 

2.4 Quadratic Expansions 

The development of the finite element as presented here, uses the 
unknown coefficients of a series expansion as generalized coordinates. 


However, an expansion in the form 

■ o F - + ±0 


(2.4-1) 


suggests the possibility of extending the expansion in q's to include 
higher orders (i.e., terms involving the square of q) . Although the 
approach adopted here does not prohibit such expansions, it does make 
them seem rather unnatural. However, such is not the case when the 
problem is approached from a physical rather than a mathematical view- 
point. It is quite natural to interpret a linear displacement field on 
a surface as a rotation of the surface. However, if rotations are taken 
as generalized coordinates, the displacement is not linear. 

All of the remarks above would be largely academic if the linearized 
equations resulting from a quadratic approximation were the same as 
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those resulting from a linear expansion. However, the essentially 
Lagrangian approach adopted here results in a different system of 
equations when quadratic terms are included. The best that can be 
hoped for is that the systems of equations resulting from the two 
methods are mathematically equivalent, although non-identical. The 
two systems are equivalent (i.e., have identical solutions) if each 
equation of one set is an independent linear combination of equations 
in the other set. If the equations are equivalent, then inclusion of 
quadratic terms is optional, and, for convenience, they may be 
omitted. That the resulting equations are equivalent, has been 
demonstrated both analytically and by example. Appendix D.l presents 
an analytical demonstration of the equivalence of the equations. 
Appendix D.2 presents a particularly challenging example of a 
gimbaled mass for which the nonlinear field is more natural, and for 
which the nonlinear system of equations Is more easily derived. It 
is also noted that the same result holds for a sequence of rotations 
about inertially fixed axes, as for the body fixed sequence in the 
example. Based on these results, it may be concluded that quadratic 
terms may be omitted. Therefore only linear terms will be retained 
in the rest of this presentation. 

2.5 Substructure Modes 


The motion of a substructure, isolated from all other substructures, 
is governed by the three scalar second order equations in R (2.3-3), 
by the three scalar first order equations in ft* (2.3-4), by the 

three scalar first order equations in C* (2.2-13), and by "N" scalar 

p 

second order equations in ^ (2.3-5). Unfortunately, "N H is frequently 
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large, sometimes hundreds or even thousands. Moreover, much of the 
motion governed by these equations is of the uninteresting, small, 
high frequency type. It is desirable to reduce the dimension of the 
problem by eliminating the uninteresting motions and retaining only 
the larger, usually lower frequency, motions. Modal analysis provides 
a basis for this reduction by uncoupling the equations and providing 
the frequencies associated with the degrees of freedom. The frequencies 
provide one of several possible criteria for truncation. The uncoupled 
nature of the equations guarantees that the solution of the remaining 
equations is unaffected by assumptions on the truncated modes, and 
allows formal establishment of bounds on the truncation error. 

Calculation of modes is possible only for constant coefficient 
linear equations. The substructure equations are inherently nonlinear, 
and even if linearization about a nominal motion is possible, the 
coefficients are not generally constant. However, the need to reduce 
the problem is sufficiently overriding that the method is frequently 
applied to problems for which it is not formally correct. In the 
presence of good engineering judgment, this has been quite successful. 
Before pursuing this tack (in 2.6), however, it will be helpful to 
identify the class of problems for which the method is formally 
correct. 

The equations will be linearized by assuming that the . q represent 
a small disturbance about the nominal motion (the rigid body motion 
solution when all ^q - 0). (It should be noted that there are cases 
for which the nominal solution is grossly different from the frame 
motion in the actual solution.) Now, arbitrarily partitioning the 
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force Into a nominal force , J?, and a disturbance force , ^F, the 
nominal motion equations corresponding to Equations (2.3-3, -4 and -5) 
are 
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Similarly, the disturbed motion equations are 
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r - li-w ^\}r ± T 

is 


A 


j 


Q . 


(2.5-6) 


Equation (2.5-3) may be written in integral form as 

,2 


k 


"-'1^7 - + j &- x o t ,+ a’x(n'x^) 
m f* • *| + 3 q • 


(2.5-7) 


Noting that in the limit as the mesh size approaches zero, the 
become a complete function space, it may be concluded that the body 


68 



force must satisfy 


X - $ R + Q f x qP» + g* x kv x qPA +°X , 

_8t * ' ' 


(2.5-8) 


where }C is constant in time and equilibrates the initial and boundary 
stresses. It may be noted that this body force also satisfies Equations 
(2.5-1 and -2). 

This leads to the perhaps surprising conclusion that linearization 
may be about an arbitrary time varying nominal motion, if the nominal 
motion is driven by body forces in the form of Equation (2.5-8). It 
may be noted that the prospects of finding a physical source for body 
forces of the indicated form are rather dim. However, this result 
suggests an analytical technique. Any forces applied to, the sub- 
structure may be represented as the sum of the equivalent body force and 
a' disturbance composed of the applied force and the negative of the 
equivalent body force. This technique provides a formally correct 
solution and reduces the structural equations to linear form (2,5-6). 

There is one potential difficulty with the method. The linear 
solution is predicated on the $ remaining small. For those problems 
in which flexibility grossly alters the nominal motion, the gross 
difference between the nominal and actual motion would appear in the 
rigid body modes inherent in the (The rigid body modes are a con- 

sequence of the unconstrained, displacement method, finite element 
formulation adopted here.) There are two ways in which this difficulty 
might be overcome. If this behavior is anticipated (e.g., by “energy 
sink" analysis) or is detected during integration, a pair of body 
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forces, a nominal force, and an equal but opposite disturbance force 
(summing to no net force on the system), could be applied to transfer 
the rigid body motion in the ^q to the reference frame variables and 
maintain the smallness of the <^. Alternatively, when the motion in 
the rigid body modes exceeds some predetermined bound during inte- 
gration, the motion in the rigid body modes could be transferred 
to the reference frame, and the rigid body motion in the _^q reset to 
zero. 

The presence of rigid body degrees of freedom in the .q, is a 

j 

consequence of the displacement method. The method also results in 
redundant Equations (2.5-4 and -5). Equation (2.5-6) is sufficient 
to determine the ^ . Qualitatively, this is a system of ordinary 
differential equations of the form 

Mq + Cq -I- Gq + Kq + Jq + Aq * Q . (2.5-9) 

The terms have the following characteristics. The "mass'' matrix, M, 
is constant, real, symmetric, and positive definite (in some special 
cases positive semi-definite). The "damping" matrix, C, is constant, 
real, symmetric, and positive semi-definite. In the special case of 
proportional damping, it is assumed that C can be expressed in a Caughey 
series. The "coriolis" matrix, G, is constant only if is, other- 
wise it is a function of time. It is real but skew- symmetric. The 
"stiffness" matrix, K, is constant, real, symmetric, positive semi- 
definite and consists of two parts, the elastic stiffness and the 
"geometric" stiffness associated with the preload. (If the preload 
depends on spin, then K is not constant unless is.) It may be 
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noted that in the ease of structural instability (buckling) K 
becomes indefinite. Such cases are excluded here. Sometimes K is 
made complex to incorporate ’’structural" damping. The "centripetal" 
matrix, J, is constant only if is, otherwise it is a function of 


time. It is real, symmetric, and negative semi-definite. The 


"spin-up" matrix, A, is null if the spin is constant, constant if 
3*o 

SV is constant, and otherwise is a function of time. (In some 
formulations A is a function of spin; see 2.5 and related appendices. 
The spin independence is a result of the Lagrangian formulation.) It 


is real but skew- symmetric. 


The equations, in general, must be integrated numerically. 
However, if the coefficients are almost constant, certain stability 
results may be derived, and if they are periodic, Floquet theory may 
be used. However, in order to calculate modes, they must be constant. 
Consequently, the nominal angular velocity must be constant. In 
addition to systems for which the spin is constant, it is customary 
engineering practice to approximate a slowly varying system by time 
slices (i.e., assuming that the deviation from constant coefficients 


is negligible over some period of time). This can be formalized by 
selecting the nominal body forces to provide constant spin, and 
absorbing the change in the disturbance equation for some period of 
time. Although a large constant spin and a small constant rate of 
change of spin might be assumed, for formally constant spin A vanishes. 

In general, second order equations must be converted to first 
order to obtain modes. Appendix A. 15 outlines the method. The 
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special case of zero spin and proportional damping is an exception for 
which eigen solutions can be obtained directly from the second order 
equations. 

When a substructure must ultimately be connected to others, a 
number of researchers have noted that "free" modes of the substructure 
do not always provide the most rapid convergence during a process of 
modal synthesis. A number of alternative procedures have been 
suggested which may speed up convergence. These include constraining 
certain degrees of freedom, adding masses and adding stiffnesses. A 


comparative analysis is forgone here; it is merely noted that all such 
methods result in a transformation of the variables, and to 
a smaller dimension set of coordinates, rj. 



(2.5-10) 


and a reduced order set of equations 


[a] (n) + [8] (o> = W . 


(2.5-11) 


The analyst's choice of method is based upon a number of factors 
including: ease of calculating the modes, storage requirements for the 
reduced equation, rate of convergence, and ease of coupling sub- 
assemblies . 

As an example of the modal analysis procedure, the component 
mode method, generally attributed to Hurty (i.e., constraining certain 
degrees of freedom) is applied to a very simple substructure in 
Appendix C.3. 
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2.6 Extended Modal Analysis 

As noted in the previous Section (2.5), the benefits of modal 
truncation are sufficiently important to justify application to 
problems for which it is not formally correct. In such cases, even 
though linearization is not formally possible, a nominal state or 
perhaps several nominal states are chosen to linearize about anyhow. 
Similarly, constant coefficients are selected, and modal analysis 
proceeds as before. The extended procedure now departs from the 
method previously outlined by relaxing all the assumptions and sub- 
stituting back into the governing equations. The Equations (2,3-3 
and -4) assume the form 


M 


9 1' 


R + £' x * qP’) 


t . T q 

+ .a n + ,3 4 n =* F 
3~ j j“ j “ 


ft r + ft 1 x S I • n» + S M ?F’ x~ 

ot ~ — — — U“ 


9t 


j r + R a n + r 3 .n = s t , 

2 “ 3“ j J” 3 

( 2 . 6 - 1 ) 


and after premultiplying by the transpose 
Equation (2.3-5) becomes 

- «'• S .I • Q'+ S M ®P’- R + S 0 
- j- “ j" 9t 2 - J- 




«'+ 


M 


t rr oi -r , . Ot . 

9t — jk k 


* + > k n - 

( 2 . 6 - 2 ) 


If the complete set of modes is retained, there is no compromise in 
the above equations. They contain exactly the same information as 
their non-modal counterparts. The a and 3 coefficients are functions 
of the variables R, , and £2* , and the equations are no easier to 
solve than their predecessors. 
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At this point engineering judgment is invoked. It is assumed 
that the effect of the modes selected for truncation on Equations 
(2.6-1) is negligible. It is assumed (as in normal modal analysis) 
that the response in the truncated modes is insignificant in the 
result. And, it is assumed that the coupling in Equation (2.6-2) 
between the truncated and retained modes has no significant effect on 
the retained modes (and conversely that the retained modes do not make 
the truncated modes important in the result). Under these assumptions, 
the modes may be truncated and the reduced set of nonlinear equations 
integrated. It should be noted that these assumptions are frequently 
reasonable and lead to good results. Some truncation can almost 
certainly be justified since the highest modes are usually more a 
reflection of the details of the finite element model than of the sub- 
structure itself. The degree of truncation suitable for a given 
problem is sometimes determined by repeating the simulation several 
times, with differing degrees of truncation, and observing the rate of 
convergence of the solution. 
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Section 3 


STRUCTURES 

A structure will be defined as a finite collection of interconnec- 
ted substructures. The interconnection mechanisms will be classified 
in three categories: those producing known forces, those producing 

known relative motions, and those with a known relationship between 
forces and motions. The interconnections producing known forces will 
be treated as external forces on the system. The interconnections pro- 
ducing known relative motions constrain not only the nominal motion of 
the structure but the disturbed motion as well. The interconnections 
yielding relationships will be- assumed to be in the' special form of the 
sum of a large known motion (i.e., constrained for nominal motion) plus 
a linear force-displacement relationship for small disturbed motion. 
Although arbitrary connections can be treated by direct numerical 
integration of the connection relations and equations of motion, a 
reduction of the equations is possible in certain cases. These reduc- 
tions are the primary subject of this section. 

3.1 Constraints 

An interconnection relationship that prescribes some aspect of the 
relative motion between two substructures is called a constraint . Con- 
straints may be broadly categorized as either equality or inequality 
constraints. Inequality constraints will not be specifically treated; 
rather it is recommended that they be treated as two separate problems, 
unconstrained, and equality constrained, with the integrator constructed 
to detect the passage from one regime to another. 
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There are some discrepancies in the literature as to whether or 
not an inequality constraint may be classified as holonomic. However, 
among the equality constraints, a constraint is classified as holonomic 
if it may be expressed in integrated form 


2 5 N S, t) = o, 


(3.1-1) 


and non-holonomic if it must be expressed in differential form 
N 


2 d 1 ? + a dt = 0, 

i=l 


(3.1-2) 


1 i 

where the £ are generalized coordinates, and the coefficients a, a 

are functions of and t. The differential form of constraint 

(Equation 3.1-2) is known as Pfaffian form . It should be noted that a 

holonomic constraint in the form of Equation (3.1-1) can always be 


written in Pfaffian form by setting , a = 


i * 


3d) 

* ; but the converse 


is not true. 

Constraints may be further classified as catastatic if a = 0 and 
acatastatic if a f 0, A constraint is scleronomic if it is independent 
of time, and rheonomic if it is explicitly dependent on time. Again 
there are discrepancies in the literature as to whether this classifica- 
tion can be applied to non-holonomic systems. 

3.2 Interconnections 

The position of any point in a substructure is defined by the 
variables of C^, R, and q , Thus a non-holonomic (Pfaffian) constraint 
between the s-th and t-th substructure may be expressed in the form 

s A s n a,s„3 , t A t_a,t_3 , s-. jS_. , t_ .t^ s r ,sa 
Cd C + A 0 C d C + B • d R + B • d R + I d,q 

ap y y ot3 y y - j 


+ C r d^q + Adt « 0 
J 3 


(3.2-1) 
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The variational form governing the structure is the sum of the 
forms for the substructures and terms governing the interconnections. 
The constraints may be incorporated via Lagrange multipliers. If there 
are "C" such constraints among the "S" substructures, then there are 
"C" subsidiary conditions of the form 

£ ( cs a • v + cs b • s r + c *r ^ ?q + C A - o ) (3 - 2_2) 

s=l ' 3 3 / 

where c - 1, •’**, C and 

,e 3 fs~a cs 


_ cs A 3 


= e 

In addition, the term 
C S 


CIB 3 s x ’y 


-/.* E E C W 

J T c=l s=l ' 


C ' ‘ CS A . S C“ <S S C 0 + CS B • 6 S R + CS 
a3 y y 


?r 6 s 5) 

3 3 ' 


(3.2-3) 


is incorporated in the variational form. This results in an additional 
term in the translational Equations (2.3-3) for each of the substruct- 
ures 


iL , i_ So, x 8 , , 11 s„, x S-, BA 


S nl SA 


,P r . q + 2 0 fl’ x “p» Tq 

j- J - j- J 

+ V x (V X »p.) + V X X 


S F + E C X CS B 


at 2 ^ at - o- T at - 

c 

5 

c=l 
(3.2-4) 

The substructure rotational Equations (2.3-4) also gain an additional 
term 
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41 V+ s a' x S I 


s’ . 9 s n ’ s T s* 

at - ' - = - at - ‘ j- j q ' j- at - j 


i T a” , S T 3 S 0 ’ SA 

I <1 + J ‘ F 2 < < J 


+ S s' • ji + ji • S s’ + V * ( v • ji) ?q + v * a i • y 


SA 

^<1 


S M S'xi S R+ ®P'X £ 


S„ SA . SA 


sA 


fl- 


at 


2 - 


j- 


3t‘ 


- 5 + “0 “q + ®0 X JO Jq 


= s t + £ 

c=X 


C X CS A 


(3.2-5) 


Finally, the deformation Equations (2.3-5) acquire an additional term 


a n a -r 

- 2. ji 


V, - 8 n\ . s 


sA 


_ So s/t , a S— W s c 

• 2 jj, q + M j- # — 2 R + 2 E • 

at 


S - SA 

y-1 1« 


sA 


s ,, sA 


, s n a T Sqi , s „ a ’ S , A S M S A s _ oo o-- 

+ jO • ■§t D + jy2 * 3t 2 * q + jJ, M * q + jJt C Jl q + JS, K !L q 


- ?Q + £ 

j c-l 


7 


( 3 . 2 - 6 ) 


In general, these equations or their counterparts in modal coordin- 
ates (i.e., with Ti’s Instead of q’s as in 2.6) with Equation (2.2-13) 
after some action to eliminate the redundant rigid body modes, must be 


numerically integrated for the X as well as the C^, R, since 

c V a B-F are > in general, functions of S Cg, S R, ^q, and time. If there 
are force deflection relationships, the forces become additional param- 
eters and the relationships become additional equations. 

The redundancy between the rigid body modes Inherent in the q’s and 
the frame variables permits the partitioning of each constraint into a 
nominal and a disturbance constraint. The nominal constraint is deter- 
mined by setting all _^q =• 0 (i.e., the zeroth order term in ^q) and the 
disturbance constraint is what remains (i.e., the first order terns 
in ^q) . Although doubling the number of constraints is not usually 
possible, the redundancy permits it here as long as the nominal 


aft* sA 
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constraints are not contradictory (e.g., no rotation between two 
substructures about an axis and a prescribed rotation about a parallel 
axis is not permitted). This is equivalent to not allowing the inter- 
connection to enforce deformation. This may be relaxed, if the 
enforced deformation is sufficiently small, by including the incom- 
patible term in the disturbance constraint instead of the nominal con- 
straint, If the _.q are thought of as the primary variables, and the 
frame variables as redundant, and the distrubance constraints are 
thought of as primary, the redundant nominal constraints can do no 
more than eliminate some of the redundant frame variables. If the 
constraints are not contradictory, any excess of nominal constraints 
can be at most redundant. 

The doubled set of constraints may be introduced via Lagrange 
multipliers, as before, except there are now twice as many, and their 
coefficients are either zeroth or first order in the _.q. The analysis 
proceeds by eliminating constraint forces, excess coordinates, and 
redundant equations wherever possible, recomputing eigenvalues and 
vectors for any constant coefficient portions of the equation, and 
integrating, taking some action on any remaining redundancy between 
frame variables and rigid body modes. The following sections outline 
the procedure in more detail for a special class of constraints. 

3.3 Rotational Constraints 

The equations governing the behavior of a structure may be sim- 
plified for special classes of interconnections. The remainder of 
this study will illustrate the procedure for one such class. Each 
interconnection will allow only one translational and one rotational 
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degree of relative motion along a single axis passing through and con- 
nected to a single point, O', on a substructure, s, and a single point, 
T, on a second substructure, t. (Any substructure may have several 
such interconnections at several points.) Thus four degrees of free- 
dom are completely constrained and the remaining degrees of freedom 
may be either unconstrained, a known function of time, or a small 
motion satisfying a known force-deflection relationship. The restric- 
tion to a translation and a rotation about a single axis (known as a 
screw displacement) is not as restrictive as it might appear, since 
trivial substructures (i.e., massless rigid bodies) may be introduced 
to represent more complex interconnections. 

Two substructures, point connected with one degree of rotational 
freedom must have bases related by 


a a 3 


‘3 3 3 

*s 1 \s 2 * »s 3 
l d y 9 y 9 y J 


*t 1’ t 2%t 3 

o y 3 y ay 


10 0 
0 cos St 0 sin 

0 -sin st 0 cos 


st. 


St 


(3.3-1) 


whereby is a dextral orthogonal basis in substructure s, and 

3 s y 

3 3 3 

—t— , in substructure t, with both and - — ■■■ directed along the 

3 y 9 y ° y 

n st " 

axis of rotation from s to t. 0 is then the rotation of t relative 

to s in the sense indicated. Such a set of axes can always be defined 
in terms of the deformed body axes by a set of direction cosines 
depending only on the substructure geometry 


' 3 sg 3 __3 

8 sd x ' 01 « a s y 6 


(3.3-2) 


> 
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(The image of the material event space axes under the position mapping 

provides such a set of coordinate axes which participate in the deforma- 

g 

tion.) The direction cosines between these axes, — - 3 — r- , and the 

9 S £1 Ut 

x 

g 

coordinate axes of the substructure, r~ , are given by 

3 s x 01 


9 sd.e 3 
T7T C. 


~s a 
3 x 


where 


"a ~sd 
3 x 


3®P’(a) a 


sd.3 P a 1- ^ • sa 

C = 60 + J ‘ 0 r .3 


(3.3-3) 


H g-e J 

when the images of the material axes coincide with the physical axes 
in the reference state (otherwise there is a constant matrix which may 

Q 

be absorbed in the geometric, . And where the direction cosines 

have this simple form because the elongations and shears have been 
assumed small in comparison to one. Finally the substructure basis is 
related to the inertial basis by 


3x 


s_$ 9 

a a a s f 3 
d x 


(3.3-4) 


Combining these relationships provides 


sg c 3 sd c Y = st c a tg c a td £ t n ^ 


y 

or rearranging 


st c a _ sg c 3 sd c Y s c 5 t 5 t:d c n tg c a 

$ Y 6 e e C n 


Expanding the deformations gives 

Y 


st a 

C 3 


sg c 3 


(s.&iSL 

\ y a s xY 


a A) s C 6 ¥ 5«+ „ 

1 7 £ e V * 3^ 


(t)^ 

_j 


(3.3-5) 


(3.3-6) 


"q ) t£ C a 
j / n 


(3.3-7) 
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The nominal constraint (i.e,, with ^q = 0) may be obtained from this 


stO a _ sg $ s y t r e tg a 
p y 0 6 e 


(3.3-8) 


This relationship may be used to determine either or S C^ in terms 
of the other and ^0. The equation is not suitable for direct inclusion 
as a constraint because the direction cosines have only three independ- 
ent terms (i.e,. Equation (3,3-8) represents three constraints) but 
Equation (3.3-8) represents nine equations. The number of equations 
may be reduced to three by writing the Pfaffian form 


1-6 st 0 r Y _JL stO r a 
~ 2 r a 3t * 


(3,3-9) 


Inserting Equation (3.3-1), it may be seen that the values of Equation 


(3.3-9) are 6 = 1, °0 ; 6 = 2, 0; 6 = 3, 0. Defining 


st0 o ’ stO* 3 

M — t) — i-^ 

a 1 ^ 1 


(3.3-10) 


where 


_J t g g a 

)V a " a 3 t0 y B 


and inserting Equation (3.3-8) in (3.3-9) gives 


st V 6 


. Ic 5 |sg r 3 s e t n tg y! 

-zScL c e c ; c ? c nJ 

Tsg 8 3 s 9 t < tg a sg 0 s r 0 _J_ t k tg al 

I l e7t S S k + c e S at c i C <J 


(3.3-11) 


Noting the identity, = 6 a \ and regrouping terns gives 
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(3.3-12) 


‘V 4 ■ (■'! K '<){- i < X '< ‘i ‘ c " X) 

* ( t£ I 5T '<)(- K X X) 

which is in the form of Equation (3.2-1). The associated vector 
coefficients (3.2-2) may be calculated as „ 
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(3.3-13) 


Noting the identities in Equations (A. 9-15 and -17) of Appendix A, 
and using Equation (3.3-8) these formulae become 
s6 . 
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(3„ 3-14) 


and 
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(3,3-15) 


The computation may be reversed 
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(3,3-16) 


since the A are clearly skew- symmetric . Inserting this expression 
in Equation (3,3-12) provides 

St % ^ = " *2 • S6 A - * t6 A (3.3-17) 

since S C^ “ S C 3 and S< ^A 3 are both the betath component of 

3 

vectors expressed in the r basis. Noting Equations (3.3-14 and -15) 

3 s x 

and rearranging 
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(3.3-18) 


^ S o f Q t 

This provides either ^ or ^ in terms of the other and 0 (much as 
Equation (3.3-8) did for direction cosines). 

Returning to Equation (3.3-7), the disturbance constraint 
obtained from the Pfaffian form analogous to Equation (3,3-9) is (to 
linear order) 
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(3.3-19) 


Since the nominal constraints must be satisfied, any linear combina- 
tion of products of the constraints with terms of linear or higher 
order in 9 q may be added to the disturbance constraint. In particular 
the second terms may be replaced by 
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Using the reverse of the identities, this may be re-expanded to 
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(3.3-21) 


This may be recognized as the negative of the middle part of the first 
term in Equation (3.3-19). Defining 


s t * 

n 


st e 3 


*t 1 

o y 


(3.3-22) 


Equation (3.3-19), as modified by Equation (3.3-20) may be written as 
st a r 6 _ stO'p tg p td it tg 5 


TT 


= _ 1 e 6 r S g c B .d e s c? t 0 td t tg y] 

2 ya L £ Cno 0 ij 


fsg c 3 *^( a )X sa s X tv td £ tg a 
L k i p p v C 


+ sg c B sd c£ s c* V tg c“l 


S<$ r SA . tdp, tA 

i r i q + / j q 


(3.3-23) 
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where 


o 


(a)' 

s , i 

i* (0) e “ 

0 X 


which is in the form of Equation (3,2-1). This may be expressed in 

9 

vector form by multiplying by ■■ and summing 

8 t 0 

y 

st n '- st0 n'P tg c p td c w tg c 6 -i-x 

- O O TT 0 

•9y . 

st | stO J 

“ '2 - 9 

- - i E 6 [ sg c e sd c E s c c t c 0 td ci tg c Y l 

2 ya L e C n n © ij 

r sg c e ^(a) A “4 s c X V td c 5 tg c“ 

[ k i Y k i M M M v K 

+ sg c B sd < S C A V t 8 c“l 4 -r 

K X p u r Cj €J 9 t y 6 

- c ?i W Ki’p E cj)^ 

,i t,, ^ tA / i ,o \ a 

«(-T £ «c)^ 


+ td c 


2 l a ii - 2 "is j nx; 5 j g |- 


- — e * .q - “ e .r ^(t) 


(3,3-24) 


The product (tyq may be recognized as ft for that direction cosine matrix. 
However, since the rotation is zero to zeroth order in q, the angular 
velocities are defined by 


1 *™° - + T e ?x - + 1 z lx (*-lL l it) 

\ 3 S ^ / 


(3.3-25) 
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which are just the theory of elasticity definitions of the rotation 
due to deformation. It may be noted that these are known constants 
once the assumed deformation fields for the finite elements are estab- 
lished (in the case where nodal rotations are selected as generalized 
coordinates, and if the points of connection between substructures are 
nodes, this is just a Boolean array, one for the i = the rotation at 
the point, and otherwise zero.) Defining the vectors 

3 x 


j£(t) - ^(t)° ~t 0 
J 3 3x° 


(3.3-26) 


Equation (3.3-24) becomes 


st. 


st*0. 


2* “ BU 0 “ “ Ccr) + ji(T)^q (3.3-27) 

The terms T in Equations (3.3-23 and 3.2-1) may now be recognized as 

s6y, S . , l 3 

i r “ " 7t ~ 


3 y 

td^ t,, , , 3 

j r = j* ( t) ^7 • 


(3.3-28) 


This disturbance constraint may be thought of as providing either some 

of the Jq or the S q in terms of the other and St fi' and 

J i — - 

In summary, a one axis rotational constraint is incorporated via 

Lagrange multipliers with constraint Equations (3.3-12 and -23). The 

additional Equations (3.3-18 and -27) are incorporated in the system 

of equations with the definition in Equations (3.3-10 and -22) implicit. 

g|*Q 

The nominal rotation, 0, is either a known function of time 
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(including zero) or unconstrained. In the latter case, the first 
line of the constraint equation may be though of as introducing the 

new variable, St ^8. In this case the variations ^A^6 st ^0, - A 6 st ^0, 

1. 1 stO stO st 

and X <5 0, must be Included. When 0 is unconstrained, 0 may 

be selected in two ways. First, if it is desirable to include the 

total (nominal and deformation) rotation in the frame state variable, 

s tn st q t n 

0, set 0 - u 0. In this case the first line of Equation (3.3-23) 

constrains the deformation rotations to be the same in each substruc- 

1 1 stO 1 1 s t 

ture. In this case, the variations - A 6 0 and A 6 0 cancel, pro- 

0 1 

viding the stationary conditon A = 0. Second, if the frame state 

s 1 0 s t 

variable, 0, only represents the nominal part of the rotation, 0, 

is arbitrary. In this case the first line defines the new variable, 

st 01 11 11 

0, and there are stationary conditions A - A and A =0. The 

second method is preferable because the associated free boundary is 

consistent with the physical system, and hence modal convergence may 

be improved. (The first method does, however provide an opportunity 

for modal synthesis). 

s 1 0 

When the nominal rotation, 0, is a known function of time 
(including zero), the rotation due to deformation may be either known 
or arbitrary. A known function may result from several situations. If 


the nominal motion is a non-zero function of time then the deformation 

s t s tO 

may be the correction for the deformed axis _0 - 0, or the sum of 

such a term and a small known motion not included in the .nominal. If 
the nominal motion is zero, the deformation may result from a small 
known prescribed motion (possibly incompatible with other nominal 
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constraints) or may be zero. Such cases can be directly included. 

However, the case of arbitrary motion requires special treatment. 

Any arbitrary disturbance motion (in the presence of a known 

nominal motion) must be small. To maintain smallness there must be, 

in general, a restoring force. If the force lies outside the system 

being treated (e.g., a controller), the disturbance rotation, 
s 1 1 s t s tO 

0- 0- 0,can be introduced via the first line of the constraint 
and it may again be noted that ^A^ = 0. The restoring force is treated 
as an external force on the system. If, on the other hand, the restor- 
ing force mechanism is included in the system, there is a known force 

s 1 1 

deflection relationship governing 0 and the Lagrange multiplier 

does not vanish. If F = f^ St ^0 ) is the force on the mechanism then 
Stl 

the term Fd 0 must be included in the system variational expression 

stl 

to reflect the virtual work done on the mechanism. If F6 9 is a 


perfect variation, the form 6V may be included instead (e.g., if 

F = K St ' l '0, V « K^ St ^0j^). In either case, the coefficient of the 

stl 11 11 

variation, 6 0, is now A + F with the stationary condition A = -F. 


The substitution 


-f( stl e) 


may be made for A in the equations, or 


^A^ may be retained as a variable and the equation ^A^ = - f ( St ^0j 

included in the system of equations* (With — - directed from s to t, 

11 9y i 

A is a clockwise torque applied to t and s is a clockwise 


def ormational rotation as viewed from s.) 


3.4 Translational Constraints 

This section continues the subject of single degree of freedom 
point constraints, introduced in the previous section. In this sec- 
tion, however, the single degree of freedom will be taken to be 
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translation along the axis 


(\S 1 r>t 1 

3 y o y 


^ . If st ? is the magnitude of 


the translation, then defining 
st^ = st^ 3 


a t 1 

3 y 


(3.4-1) 


The position of the connection point of t may be expressed in two 
ways 

C R + t "P (t) = S R + S P (q) + St g (3.4-2) 

or rearranging 
st 


I = + qF ? (t) + jP(t) Jq “ S R » qP ( 0 ) - ®P (a) ®q, (3.4- 


3) 


The nominal part of this is 


st °5 - h + ^' (T > - S E - o£-' (a) 


(3.4-4) 


or rearranging 


t R = St0 i + S R + qP(<J) - qp'(t) (3.4-5) 

t s 

This relationship may be used to determine either R or R in 

s 1 0 

terms of the other and £. This is not in Pfaffian form. The 
Pfaffian form will now be obtained for inclusion in the system equa- 
tions. 

9 

Taking the dot product of Equation (3*4-3) with — - — gives 

a t oc 

y 

Stf m tg c a td c B t c Y t R 6 + tg c a td c B t p > (T)Y 

St a sg 3 sd Y s 6 s e _ st a sg 0 sd y s ’ 


, Q w C’ C R 

3 y 6 e 


B 6 c; -p (a) 


- tg c“ td c* ^ V + tg c“ td C® V (T)^ 

_ tg a td J t Y s 6 tg a tdj t y s e s ’ .i 

3 y 6 R y C 6 C 6 P {0) 

where Equation (3.3-6) has been used. 


(3.4-6) 
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Differentiating this gives 


.. E « . t*. tj y , t, c „ .» y (X)Y 

- “eg ,w s ;s v e •«* - '*<» ‘tmi “S ‘cj -c« v«>« 


'3 r v 7 3 j 


tg a td 3 _3 t if tg a td p 3 t Y s 6 

+ C 3 3t <5 k " l 3 Y 3t 6 K 

tg a td 3 3_ t y s e Sp 1 / n £ tg a td 3 t y 3_s e Vfrrt 6 

C 3 C y Jt C <$ C 6 P w C 3 C y C 5 3t C 6 P {0) 

t g _ot td 3 t y _3 t 6 tg a td 3 t Y _3 s <5 

+ C 3 Y <5 3t R " C 3 Y 6 3t 

+ *«<* td cJ Jp’(t)Y ^ - tR c“ td c® ‘cj s Cg ®p’ (a) e . 

(3.4-7) 

The nominal constraint may be determined from this, or by differentia- 

9 

ting the dot product of Equation (3*4-4) with ■■ ■ ,■— , as 

9 to c* 

y 

a St0 r a tg a 3 t 3 t <5 t i t e tg r a 3 t p 3 t 6 t 6 s_e 

3t” 5 C 3 3t~ C Y C Y C e R ” C 3 "at" C y C Y e R 

t 8c a 5 t c.^ t r (S s n^ s p\nY > - tg n a 3 c 5 s n^ -£ s n^ 


'6 3t 


« - c w - C » -p (a) 

Y Y e e 0 v 


’C~ d C u ■£ nP 

3 y Y e ,9t e u 


tg a t B 3 t v tg a t B _3_ s y 
+ C 6 y 3t R C B Y 3t R 

- V si- M nr 1 '!) 
W»r V )*\(j?- V ) 


(3.4-8) 


This is in the form of Equation (3.2-1) the associated vector co- 
efficients (3.2-2) may be calculated as 
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sa A = sa A 9 


6C f n 

ox 


- / tg c“ ‘c® S c s — ,.V ft'w' — r r 

I 6 Y Y 3 s x '<5) 1 0 3 s x ' ? 


x > (a) 


.to *a 0- 

3 y 


and 


(3.4-9) 


ta. n<53 ta. 3 

A = e A 


5 ®^ 


E n66 tg c a t c 6 t R e _ tg c « t( ,6 s R e 
Y p £ p e 


tg a t 6 s j; s ' ( xC\ 3 

C 3 e e 0 P {a) ) 9 t x f n 


t8c ®^’® ) x \< tR V 

t d x / \ 3 x 


J! \ + ( tg a 3 \ x /t J s e S_ 

VV \ ® 3 V®/ l 6 aV 


- JL- x t R H - — r X S R + 5__ -X S P (a) 

.to ? a - .to a - .to T a 0- v ; 
o y 3 y 3 y 


(3.4-10) 


Similarly, the other coefficients are 


B - - 


tg c a t c 3 3 

6 Y 3x T 


.to 'a 

3 y 


and 


'B •» 


tg c“ V -i. 

® Y 3x Y 


.to *a 

3 y 


(3.4-11) 


(3.4-12) 
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which may be recognized as the nominal part of Equation (3.4-2), dif- 
ferentiated in the - • - ■ — reference frame. This may be rearranged 

d tu 

y 

slightly to 

C R = ^ x fc R .-‘a xVV ®p' (a) + S R + s n x Sp'(a) + st % 

(3.4-15) 

Noting Equations (3,3-18) and (3.4-5), this may also be written as 
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K 


fc R = l a Q + st V) x ( st0 £ + s r + g]?'(a)- qp’(t)- s R- qP.' (a) ) 


+ S R + s n x sp’<a)+ st0 | 


= s S2 X st0 £ - x q/(t) + s n X gp'(a) + S R + st0 Q X st0 | 


st V x Jp’(t) st0 | 


(3.4-16) 


t • S • 

This may be thought of as providing either R or R in terms of the 
other and . 

Returning to Equation (3,4-7), the disturbance constraint may be 
obtained by first multiplying Equation (3.4-8) by tc ^C^ 


tg r i td 6 tg a stQ a 

c 0 c n n ' % 


- t8 c' td c e J- V 1 V fc c 6 V 

0 r) dt Y Y e 


tg c x td c e t n t c « t c 6 s R e 


n St 


Y Y 


tg c i td c ® j- t c n V fc c 6 s c ? ®p’(a)' 


'n St 


Y Y 


£ O' 


ts cl td c e V S C 6 S C 5 Jr S C ? ®p’(a) ? 
0 n y Y e 3t e 0 


tg 1 td 0 t n t Y tg i td 0 t n J>_ S Y 
+ c e c n C Y 9t R " c 0 n c y 3t R 


(3.4-17) 


The first, second, fifth and sixth terms of (3.4-17) are the same as the 
fifth, sixth, ninth and tenth terms of Equation (3.4-7), so subtracting 
(3.4-17) from (3.4-7) provides the disturbance constraint. 
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at 


3 st c a tg a tdj tg J _9 stOr.6 

T S'" '“'Q K, K, *}«- ^ 


’3 w y Y 3t 

t8 c“ Vt )} ]k ^ ^ + tg c“ >(T)J V<T) 


Y 


- ,g «; ;*«j ;« ■*' - e s ;i •<$ ^ ■>’<•> 

tg a td j 8 t y t e t e s n s ' *n sa 

" C S C Y W C 6 C 6 c ; C c i p (a) i q 

tg a td 6 U s e s e _3_ s n s .n s a A 

- Cg c y c 6 c 5 c ? 3t C ? ± P to) A q 

+ tg c“ td c® *p’(t)Y J$ - tg c“ td c® ^ S C^ Jp'(af ^ 

‘.n ( *'c >1 ) - “‘«v ( tc S fe e =i ) ' 

? r (i«K“ r G t «) 


sa 


sa 


(3.4-18) 


This is in the form of Equation (3.2-1). The associated vector 
coefficients (3.2-2) may be calculated as 


sa 


'A * 


0eti sa. 9 

E A en TT^e 

9 x 


_ Sen 


* e 


■-( 


(-«cj tJ c; 'c ; -c| •/<„" ;;) 


a x 


tg a td j t y s _e 9 


8 V Y <5 6 a s 'e 

9 x 


c“ -“cr ~c; ~c 

■ TT - ^ - * i ? - (a) i q 

3 y 


W w ” 


(3.4-19) 


and 
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ta A = e 9£Y ta A - 


ey a t ’0 
' 9 x 


_ ^0£y / tg n a td„0 t„e S^n s„' , .n sa\ 9 

' e r c s S c c c s i p (a) i q Tr T e 

0 X 

_ / tg c a td 3 9 \ / t e s n s * vn sa 9 \ 

( c e c y??r) x { c ? c u p (ff) i q 77?) 

= ^ x i p ’^i q . 

9 y 

Similarly the other coefficients are 

sa r = _ tg c a t d c 8 t c Y s c^pV) E 


(3.4-20) 


,t 'a i- 


p (a) 


(3.4-21) 


and 


r r ■ tgc a >^)b tp, ^ Y 

- t8c g M s « 6 - tg c^ J*(T)> V s *cl V(V 
+ tg c^ td c; Jp'd)^ 

" ji (T) x - t p (x) + S R + s p’(a)|+ jP'(x) 


t ’a 


(3,4-22) 


This can be expressed in vectorial form by multiplying by — - 

a t oc 

y 


st | _ st0 | « _ j^(T)jq X S R - ^(T) x Jjp (T) + ^(t) *jq x S R 

+ j«t) ‘q x ®p'(o) + *2 x ®p'(0)®q - \i,x ®p'( 0 ) ®q 
+ j-' (T) j q _ i-' <a) i q (3.4.23) 
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s v 

where it has been noted that ^(o)^ is skew-symmetric for small 

deformations. This result may also be obtained by differentiating 

Equation (3.4-3) in the | reference frame. 

9 y 


3 ty ’ st 
3t 


§ - t | - C S x fc R - Jj(T) Jq X ** _ ^ (T ) ^ X *p’ (T) 

+ Jp’ (T)^ - S R + ^ x S R + JiCt) x S R - x ®p' (a) 
+ t n X ®p'(o) + jid) *q X ®p'(o) - ®p' (a) *q - s n 

X ®p'(a) ®q + *fl X ®p'(a) ®q (3.4-24) 

and subtracting Equation (3.4-14), Noting Equations (3.3-18, -27 and 
3.4-5) , Equation (3.4-23) may be written as 

St !.- St0 ! - - JS!<t> Jq x ( st0 C + S R + ®p ' (o)-qp' Jqx^p’ ( T ) 

+ j4(t) jq X S R + ^(t) ^q X ®p'(a) + ( S q + St0 fi') 

x ®p'(0) ®q - S 2 X ®p’ ( a) *q + jp’(-r) ‘q - ®p' (a) ®q 
- f$<T) ^ x st0 5 + st 0 n'x ®p' (a) + V(T) J* 

■J J 


S tj ' / SA 

- il < a > i.q 


j £VW j 

(3.4-25) 

The terms in this expression may be identified by noting that the 
identity 

St S = St °5 - ®p'( 0) ®q + jP'(T)^q (3.4-26) 


differentiated in the 


~t 
d y 


reference frame gives 
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St | - ( St °| - ji(T) jq x 3t0 ?)- (®p'( 0) 3 q - St V x ®jp' (a) 


+ 1 *™ ^ 


(3.4-27) 


This disturbance constraint may be thought of as providing some of 
the ^ in terms of the and St £. 

In summary, a single axis translational constraint is incorporated 
via Lagrange multipliers with constraint Equations (3.4-8 and -18). 

The additional Equations (3.4—5 and -26 or —25) are incorporated in 

g H Q 

the system of equations. The nominal translation, £ is either a 

known function of time (including zero) or unconstrained. In the 

latter case, the first line of the constraint equation may be thought 

s 1 0 

of as introducing the new variable £. In such cases, the variations 

^A 6 s £,-V<5 St °€ and ^A^6 St £ must also be included In the variational 

form. When St ^£ is unconstrained, there are two ways to select st £. 

First, if it is desirable to include the total (nominal and deforma- 

tional) translation in the frame state variable set St £ - stO^ 

In this case, the first line of Equation (3.4-18) constrains the deforma 

tional translations to be the same in each substructure. The varia- 
1 1 stO 1 1 st 

tions - A 6 £ and A 6 £ cancel giving the stationary condition 

0. 1 — s 1 0 

A =0. The second option is to let the frame state variable, 


represent only the nominal part of the translation by making st £ 


arbitrary. In this case the first line defines the new variable, 8 
There are stationary conditions ^A^ - 0, and ^A^“ = ^A^ so that 
^A^ = 0 also. The second method is preferable because the associated 
free boundary is consistent with the physical system, and hence modal 
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convergence may be improved. (The first method, however, does provide 

an opportunity for modal synthesis.) 

st 0 

When the nominal translation, ^ is a known function of time 
(including zero) , the translation due to deformation may be either known 

4 

or arbitrary. A known function may result from several situations. If 

the nominal motion is a non-zero function of time, then the deformation 

st s 1 0 

may be the correction for the deformed axis, £ - §, or the sum of 

such a term and a small known motion not included in the nominal. If 
the nominal motion is zero, the deformation may either result from a 
small known prescribed motion (possibly incompatible with other nominal 
constraints) or be zero. Such cases can be directly included. How- 
ever, the case of arbitrary motion requires special treatment. 

Any arbitrary disturbance motion (In the presence of a known 
nominal motion) must be small. To maintain smallness, there must, in 
general, be a restoring force. If the force lies outside the system 
being treated (e.g., a controller), the disturbance translation, 
stl^ __ stj. _ can k e introduced via the first line of the cons- 

traint; and it may again be noted that ^X^ - 0. The restoring force 
would be treated as an external force on the system. If, on the 

other hand, the restoring force mechanism is included in the system, 

s 1 1 

there is a known force deflection relationship governing £ and the 
Lagrange multiplier does not vanish. If F = is the force on 

S 1 1 ' 

the mechanism, then the term F6 £ must be included in the system 
variational expression to reflect the virtual work done in the mechan- 
ism. If F6 £ is a perfect variation, the form <5V may be included 
instead (e.g., if F = K St ^£, V = y K ) . In either case the 


100 



s ti X 11 

coefficient of the variation, 6 becomes X + F with the station- 


1,1 


ary condition X = -F. The substitution -f( st may be made for X 

1 1 

in the equations, or X may be retained as a variable, and the equa- 

o 11 istiX \ 9 

tlon X = - f y included in the system of equations. (With — ~ 

dy 

directed from s to t, X is a positive force on t, and is 

a positive deformation displacement of t away from s.) 


3.5 Elimination of Constraint Equations and Forces 

Thus far, constraints between substructures have been incorporated 
via Lagrange multipliers. They may be interpreted as forces ( and 
torques) on the substructures, due to the constraints. The number of 
system equations has been augmented by the C constraint equations, and 
the number of variables has been augmented by the C Lagrange multi- 
pliers. In some cases, it is desirable to retain the X as explicit 
variables. This practice may even be computationally more efficient. 
However, if the integration forces are not of interest, there is fre- 
quently a computational advantage in reducing the dimension of the 
system equations by eliminating them (and the associated equations). 
Even if the forces are desired, it may still be advantageous to elimin- 
ate them and subsequently recover them. 

The nominal translational constraint (3.4-2) results in constraint 
forces (see Equation (3.2-4)) on the s substructure (3.4-11) 


0,1 9 0,2 3 0,3 9 

" A 'TcT T T " A 'to T 2 " A 770^3 


y 


(3.5-1) 


9 y 


9 y 


and on the t substructure (3.4-12) 

9 0,2 a 

+ A 


+ °x 1 


a to T 1 
9 y 


,to T 2 
a y 


^ + °x 3 


ato T 3 

9 y 


(3.5-2) 
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and, in addition in constraint torques (see Equation (3.2-5)) on the s 
substructure (3.4-9). 


- ~X 


o, 1 a 


0, 2 3 


0, 3 3 


_^r X S p (a) - “A~ — x =p (o) - X 0 P (a) 

3 y 3 y 


3 i- v 

y 


(3.5-3) 


and on the t substructure (3.4-10) 


/o,i s , 

\ A .tO '1 + A .to ’2 
x 3 y a y 


0,2 3 ,0,3 3 

r?r + A 


.t0 '3 
3 y 


y / S R + 0 s p’ (a) - s) 


(3.5-4) 


The torques can be eliminated from the rotational equations by sub- 
tracting the cross product of the s translational equation with 

s 1 

P (a) from the s rotational equation. And, similarly, subtracting 

S Q * t 

the cross product of the t translational equation with R + ^P_ (a)- R 
from the t rotational equation. The resulting equations replace the 
rotational equations, but the translational equations are retained 
for no net change in the number of equations. (This effectively 
translates the origins to the nominal attachment point on the s 
substructure . ) 

The constraint forces can be eliminated from the translational 

equations by summing the equations for the two substructures. The sum 

re pi aces the individual equations resulting in a reduction by three 

* * s 1 0 

equations (and three X s). If ^ is arbitrary, the constraint has 
s 1 0 

added E, to the system variables, so another equation is needed. 

The required equation may be obtained as the dot product of the t 

0 

substructure translational equation with — - - • This contains no 

oi 3 y 

constraint force since only X is involved, and it Is zero for this 
case. 
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The nominal rotational constraint (3.3-12) results in constraint 
torques (see Equation (3.2-5)) on the s substructure (Equation (3.3-14)) 


“a 1 


a t0 ’1 

9 y 


+ °A 2 


.tO 1 2 

3 y 


+ °A 3 


.to 1 3 

3 y 


(3.5-5) 


and on the t substructure (3.3-15) 


- V 


.to '1 

d y 


- °A 2 


,t0 '2 
d y 


- °A 3 


.tO ’3 
a y 


(3.5-6) 


Summing the rotational equations for the s and t substructures results 

in the cancelling of the constraint torques. The sum replaces the 

individual equations in the system equations resulting in a reduction 

1 stO 

by three equations (and three As). If 6 is arbitrary, the con- 
st 0 

straint has added Q to the system variables, so another equation is 

needed. The required equation may be obtained as the dot product of 

3 

the rotation equation of the t substructure with — ^ — yrr. This equation 

0 1' 3 y 

contains no constraint torque since only A is involved and it is zero 

for this case. It may be noted that since the rotational equations for 

s and t are in different bases, constant coefficients for either may 

not remain so after transforming to a common basis. In such cases, it 

may not be computationally advantageous to sacrifice constant coeffici- 

sto 

ents for reduced dimension. • However, if 8 is constant in time, or 
if one of the substructures is symmetric about the - — — — j— axis, con- 

3 uU JL 

y 

stant coefficients can be preserved. 

The disturbance translational constraint (3.4-18), in addition to 
forces on the disturbance equations, results in constraint torques 
(see Equation (3.2-5)) on the s substructure (3.4-19) 
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and on the t substructure (3.4-20) 


1,1 a , 1,2 a . 1,3 a sj, , sa 

A -T-’ f I + A 7F^ + X i- (o) i q 

a y ay a y / 


(3.5-8) 


1 1 

Since the A are forces on the def ormational equations, they are first 

order in q. Consequently, the constraint torques (3,5-7 and -8) are 

second order in q and may be dropped from the rotational equations. 

The elimination of constraint forces for the deformation equations 

will be postponed to the following section on modal coupling. The 

remainder of this section will be directed to eliminating the constraint 

equations. The constraint equations may be eliminated by substituting 

the right hand side of Equation (3.3-18) for each appearance of ^2, 

and the right hand side of Equation (3.4-5) for each appearance of R. 

In addition, the direction cosine relation (2.2-13) for the t substruc- 

s 1 0 * 

ture may be deleted in favor of the implicit definition of Q 

(3.3-10) and the fact that resolving the rotational equations in a 

s 1 0 

common basis introduces 0, and makes the equation second order. 


3.6 Modal Synthesis 

The partitioning of constraints into nominal and disturbance parts, 
leads to constraints between disturbance variables, q, in addition to 
those between frame variables, discussed in the previous section. As 
in the previous section, it may, or may not, be computationally 
advantageous to eliminate the constraint forces. In some cases, the 
structure may be divided into substructures, primarily to obtain 
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smaller, more tractable models. The constraints, in this case, 
merely reassemble the pieces, eliminating all six relative motions. 

In such cases, it is generally computationally advantageous to elimin- 
ate the constraint forces. If the dimension of the substructure equa- 
tions has been reduced by modal truncation, the procedure is known as 
modal coupling or modal synthesis. 

Since the single degree of freedom point connected constraints 
are holonomic, it is more convenient to deal with their integrated 
form. Noting Equations (3.3-7 and -8) 


st a _ sg J 
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= st 0 c“ + ^(a)^ £ st 0 c“+ st 0 C® ^(x)® tA ( 3 . 6 - 1) 
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ing Equation (3.3-1), this may be written as nine explicit equations 


10 0 

0 c0 s0 
L0 — s0 c0. 


"1 

0 

0 


'1 

0 

0 

c°0 

s°0 

+ 

0 

c°0 

_0 

-s°e 

c°0_ 


.0 

-s°G 

0 

+V 

t .21 
-ip 


~1 

0 

t.3 

0 

■+■ ip 


0 

c°e 

,t 2 
+ ip 

t ,i 

0 


0 

Oft 
-s 0 

i. 




- 

-V+c° 

e°0 

s , 3 
t p -s 

°eV- 

t ,3 

- 

> 

c°9+s° 


o -V y 

s tp 3 0 - s tp 1 

V V 0 J 


(3.6-2) 


-s°0 S ^ 3 -c°0 S ^ 2 + t ^ 2 f -s°0+c°0 S ||; 1 -c O eS|j 1 } c O 0+s O 0 S i^ 1 -s°e t: ^ 1 


105 



where the abbreviations are defined as 
c0 = cos^ St o), s0 = sin| st 0^, 


c°0 = cos^ St0 e), s°0 = sin 


n> 


i.a 1 a3 sg , ; vY sa t,a 1 a$ tg , tA 

tp = 2 e Y i q > ^ = 2 e y j’K' l )g J q . 


Noting the' small angle approximations, cos E « 1 , sin E=^E, and the 
sum of angles formulae sin(xfy) = sin x cos y + cos x sin y, and 
cos(x+y) = cos x cos y - sin x sin y, the (2,2) (2,3), (3,2) and (3,3) 
elements all yield 


st n StO„ Sg , l v SA tg. *, vtA 

0 - 0 = - ®iKa) .q + ?^(l).q 

i i 3 1 


(3.6-3) 


where 
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n n 


1 n I e g~a 

= 2 % c n 
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T1 r 

The (1,2) and (1,3) equations may seem to be the linear combinations 
-cos^ St ^0^( 2, 1)J + sin^ St ^0^(3 ,l)j and 

-sin( St °e) £(2,1)J- cos^ St °0^(3,l)J of the (2,1) and (3,1) 
equations. Selecting these as the two independent equations provides 

C ®1(<(T ) 3 = cos( st 0 e) s |iKa ) 3 - sin( st 0 e) s |iKa ) 2 

*|*<t) 2 = cas( sto 0) s f*(0) 2 + sin ( st0 9 ) s ^(a) 3 ^ (3.6-4) 

The translational constraint may be obtained by dot multiplying 


Equation (3.4-26) by 


~t 'a 
d y 
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(3.6-6) 


where the abbreviations are defined by 

snt sg * , vd sa . t a tg * , tA . _ . 

P = (a) and P = xT (l) ^q. In explicit form, the 

constraints are 

st (, - St0 £, - -^f'(o) 1 S 5 + V(T) 1 tA q 

i J 3 

(V (t) 2 _ stO^ . [ co8 (»fO 0 |8g p ' (o) 2 + sin (^0 Q Jsg p ' (o) 3jsA 

( t j-p (T) 3 + st °5 c ?iHT) 2 jjq = ^cos( st0 e) s |p (a) 3 - sin^ st0 e^ s |p (°) 2 J^q 

(3.6-7) 
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The constraints (3,6-3, -4 and -7) are complicated by the 
presence of the nonlinear nominal relative translation and rotation 
terms. In general, the nonlinear constraint might not be eliminated, 
or might be eliminated each integration step. If the nominal rela- 
tive translation and rotation are known functions of time, a non- 
autonomous solution to the constraints might be found. In either case 
the calculation can be simplified by selecting the rotation and trans- 
lation at the connection point as explicit degrees of freedom (e.g., 
boundary points in a Hurty component mode analysis). One possible 
simplification is to apply a stronger constraint requiring the de- 
formantional rotations and translations (normal to the axis) to be 
zero in one substructure. This would be formally acceptable to the 
extent that it eliminated still present redundant rigid body modes. 
However, it would generally be expected that there would be insuffici- 
ent redundant variables to allow this for every constraint. This 
might be an acceptable approximation if the constraint point was not 
expected to participate significantly in the deformation. If the con- 
nection point is constrained on one substructure the constraint equa- 
tions yield a similar constraint for the other substructure. 

The greatest simplification occurs when there is no nominal rela- 
tive translation or rotation. In this case the nonlinearities are 
gone, and explicit solution is possible. It may be noted that the 
constraints in Equation (3.6-3) and the first of Equations (3.6-7) are 
always linear (although inhomogenous) as noted earlier in this section,, 
The case of no nominal relative motion frequently occurs in modal 
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coupling. The remainder of this section will be dedicated to the 
treatment of this special case. 

Without loss of generality, it may 'be assumed that StO 0 and 

stO r 

s sre zero (either by changing the geometrical direction cosines 
or linearly extrapolating the displacement field). In this case the 
constraints assume the simpler form. 

at e - -><o> 1 1 $ + ^(t) 1 

0 = - s |i|i(o) 2 + ^(x) 2 

o - - s ®^<a) 3 ®q + ^(T) 3 *q 

St C - -^'(a) 1 + ^p'd) 1 5$ 

0 = - S 8p'(a) 2 ®q + C ®p'(T) 2 *q 

0 - - S ®p'(CT) 3 ®q + C |p'(t) 3 jq (3.6-8) 

The union of the sets of ^q and may be thought of as a vector 
space. Similarly the relative motions on the left of Equations 
(3.6-8) may also be thought of as a vector space. The Equations 
(3.6-8) may then be recognized as a linear function. This may be 
written in abbreviated form as 

M = [T] {q} (3.6-9) 

Since the number of q is generally much greater than the number of s, 
and prohibiting contradictory constraints, the linear function, T, 
generally has a non-trivial null space (and the matrix has a rank 
smaller than its larger dimension) . The q can be represented in a 
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new basis consisting of a basis for the null space and a set of 
vectors whose images under T span the space of s. 

{q} = [N] {£} + [B] (C) (3.6-10) 

where 

[T j [N ] = [0] 

and 

[T][B] = [U] 

B is sometimes referred to as a pseudo-inverse. (Redundant con- 
straints could lead to zeros in the identity matrix, however, this 
is of no concern here.) Substituting Equation (3.6-10 in -9) yields 

(C) = {s} (3.6-11) 

It may be noted that a "B M basis vector added to any linear com- 
bination of "N" basis vectors is still a "B" basis vector. The dis- 
turbance equations have the generic form 

[D] {q> = {Q> + [T]M (3.6-12) 

where equations for both the s and t substructures have been 

included (and D is thought of as an operator including differentiation 

T T 

in time) . Forming the linear combinations N and B , and noting 
Equation (3.6-10) 

[N] T [D][N]{£} + [N] T [D] {bHc} =[N] T {q} + [n] t mU) 

[b] t [d][n]{£;} +[b] t [d][b](c} - Ib] t {q} + [b] t [rHx) (3.6-13) 

The term N Y A may be recognized as the virtual work done by the con- 
straint force on a virtual displacement compatible with the constraint, 

and is therefore zero. Thus the constraint forces do not appear in 

T 

the first equations. This is not true of the term B TX. In fact, 
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the second equations may be thought of as determining the A. If the 
constraint forces are of no interest, the second set of equations may 
be dropped, or if they are of interest, retained for subsequent re- 
covery of the constraint forces. Only the first equations are re- 
tained in the system equations. Noting Equation (3.6-11), they may 
now be written 

[NJ T [D][N]{£} - [N] T {q - [D] [B] {s}} (3.6-14) 

Since these equations already reflect the constraint (3.6-9), it may 
be dropped from the system equations and if desired. Equation (3,6-10) 
can be used for subsequent data recovery. Thus for each such con- 
straint, there is typically a reduction by twelve equations, six A f s, 
and six q’s (in favor of ^s). 

The use of one fully constrained connection above is illustrative, 
not restrictive. For instance, three substructures with four such 
connections (as long as they are not contradictory) can be treated in 
exactly the same manner by merely increasing the dimensions of the 
indicated matrices. And, even if the connections are not fully con- 
strained the two constraints, Equation (3.6-3) and the first of Equa- 

s t 

tions (3.6-7), can always be treated in this manner. If either £ 
s t 

or 0 is arbitrary, an additional equation is required. This is 
obtained by adding the appropriate (first or third) member of the 
second set of Equations (3.6-13) to the system equations as previously 
noted, the associated A is zero for this case, so the equations has no 
constraint forces. 


Ill 



At this point the methods outlined in Section 2.5 can be applied 
to the coupled substructures. So if constant coefficients can be 
assumed, inodes can be generated for the coupled system. If con- 
strained modes are to be calculated for further coupling (e 0 g.j 
Hurty f s method) it may be desirable to defer the calculation of the 
null space until after the static shapes have been calculated, since 
the calculation of static shapes is simplified by having the boundary 
nodes as explicit degrees of freedom, and transforming to the null 
space could result in variables involving linear combinations of 
the boundary degrees of freedom. 
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Section 4 


SUMMARY 

The purpose of this last section is to summarize the concepts and 

o 

procedures outlined in the previous sections . As a vehicle for pre- 
senting this summary, the procedure for analyzing a dynamical system 
will be outlined. 

4.1 Isolation of Structural Elements 

The first step in the analysis of any dynamic system is to define 
the system under consideration. The process involves two procedures: 
defining the boundary, and defining the interaction across the bound- 
ary. For instance, for a train on a track, the track could be treated 
as external to the system, providing a rolling constraint; or in. the 
system, with the road bed external to the system providing an elastic 
foundation. 

The second step is the division of the system into dynamical 
subsystems of separate analytical types. For instance, in the train 
example, some of the dynamical subsystems might be: the hydraulic 

coupling system, the pneumatic brake system, the turbulent boundary 
layer, the biological engineer, the electrical power system, in addi- 
tion to the mechanical structure of the train itself. Each subsystem 
might be separately analyzed, and in some cases equations governing 
the behavior of the subsystem might be obtained. In some cases a 
trivial solution might be selected (e.g. constant temperature for a 
thermodynamic subsystem if it is not particularly relevant to the 
system), or the subsystem might be excluded from the system. (This 
latter case really falls under the first step, defining the system.) 
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As in the first step, the second step involves two procedures i 
isolating the subsystem, and defining the interaction. For instance, 
the acoustic subsystem loads the structure, and the structure displaces 
the gas. The interactions now include not only interactions with 
things outside the system, but also with different subsystems within 
the system. 

This study treats a subsystem identified as a structure. The 
primary interactions between structure and other subsystems are dis- 
placement and force (including body and surface forces). Although other 
interactions are possible (e.g. thermal and chemical), they will not 
be explicitly treated. It will now be assumed that all of the dynam- 
ical subsystems except the structure have been treated, and that their 
interactions with the structure result in either displacements or 
forces. Although the discussion in Section 1 has applicability to 
general material systems, the primary thrust of this study is solid 
materials. Therefore fluids and gases are included in the dynamical 
subsystems assumed to be already analyzed; not in the structural 
subsystem. 

4.2 Selection of Substructures 

A structure may be divided into a number of substructures for 
several reasons. Some subdivision may have already been made by the 
removal of other dynamical subsystems (i.e. if the dynamical subsystem 
joined the two substructures). If there are any joints in the struc- 
ture which allow large relative motion, a subdivision should be made 
across the joint so that no substructure includes large relative 
motions internally. It may be desirable to select substructures to 
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isolate portions of the structure of different types (e.g. to isolate 
a visco-elastic material, or an orthotropic shell). The structure may 
also be subdivided to distribute the burden of analysis among individ- 
uals, or in some cases companies. Finally, and frequently most impor- 
tantly, the structure may be subdivided to realize the substantial 
reductions in computational costs available through modal synthesis 
techniques. Since the savings is dependent on the amount of inter- 
face a substructure has with others, when there is no other reason 
for picking a given division, it should he selected to minimize the 
amount of interface (e.g. a long slender structure should be divided 
into two substructures by a cut perpendicular to the long axis rather 
than parallel to it). 

4.3 Substructure Dynamics 

The next step after dividing the structure into substructures 
is to obtain the equations governing the behavior of each substructure. 
Although there are numerous ways of obtaining these equations (e.g. 
finite difference, Closed form, etc.), one of the most popular and 


powerful methods is the displacement formulation of the finite element. 
This method is developed in detail in Sections 2.1 - 2.3. Whatever 
method is used, the result should be a set of equations governing the 
behavior of the substructure in the form of Equations (2.2-13), (2.3-3 


-4, and -5) which are reproduced here for reference 
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Along with the governing equations, the equation for recovering 
displacements on the interfaces (e.g. Equation (2.4-1) must be available 

r = <£’ + i*’ ±3 (4.3-5) 

The class of substructures or dynamical subsystems treated may now 
be expanded to include any resulting in equations of the above form. 

4.4 Substructure Modes 

In many cases, substantial reductions in the dimensions can be 
obtained by transforming from the generalized coordinates to modal 
coordinates. The circumstances under which this is permissible and 
the procedure for obtaining modes are outlined in Sections 2.5 and 2.6. 
Because of the substantial savings resulting from this procedure, it 
should be carried out for all substructures for which it is possible. 

The result of the transformation is to reduce the number of equa- 
tions and to modify the numerical values in the governing equations. 
However, the basic form is unaltered. (At least when rewritten in 
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first order differential equation format.) It will therefore be 
assumed that modal analysis has been performed for all substructures 
for which it is possible, the governing equations now involve modal 
variables, and the number of elastic equations has been appropriately 
reduced. Conversely, if modal analysis is not possible for the sub- 
structure, the equations remain unaltered. 

If the interconnection points are retained as degrees of freedom 
(e.g. Hurty’s method), then no back transformation is required to 
apply the constraints as developed. If they are not retained, and the 
substructure was treated as first order, the constraints must be modi- 
fied to reflect the contact transformation. Xt will be assumed that 
the displacements at the interaction boundary are available as a linear 
combinations of generalized displacements (not velocities). 

4.5 Substructure Coupling 

At this point, there are at least two distinct paths which might 
be followed. It has been assumed, thus far, that the def ormational 
coordinates contain all the rigid body modes, so there is a redundancy 
between them and the frame variables. The first path eliminates the 
redundancy at this point in the development by constraining the defor- 
mation and frame variables. For instance, a mean motion frame, or a 
principal axis frame, or a frame attached to the material at some point 
of the substructure. If this path is chosen, a constraint requiring 
that the motion of some degree of freedom be the same on two substruc- 
tures, would be reflected by a single equation involving both frame 
and def ormational variables. 
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The second path defers the elimination of the redundancy between 
deforma tional and frame variables. Consequently the physical inter- 
pretation of the frame motion (e.g. as mean motion) may be lost. Since 
each constraint can be thought of as constraining frame as well as 
deformational motion, each constraint is partitioned into a nominal 
(zero order in deformation) and deformational (first order in deforma- 
tions) part. Each part is treated as an independent constraint. Focus- 
ing on the nominal constraints for a moment, it is clear that for sub- 
structures with numerous, redundant, connections, the number of nominal 
constraints could easily outnumber the frame variables. As long as the 
constraints are prohibited from being contradictory, the excess con- 
straints are at most redundant. This procedure leaves as many frame 
variables as there are rigid body degrees of freedom for the structure. 

In this study the second path has been selected because it tends 
to minimize the coupling between frame and deformational variables. 

Since this method tends to place motion in the deformational variables 
that the first path would associate with the frame, the chances of 
obtaining deformations large enough to compromise the linearization are 
increased. In such cases the model should allow relative frame motion 
at some joint to reduce the deformations. 

If the first path is elected, the appropriate constraints may be 
written, (including those defining the frame) and the resulting equa- 
tions integrated. The following sections, however, are based on select- 
ing the second path. 
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4.6 Point Connections 

As noted in Section 3.6, rigid connections provide six linear 
constraints, and single axis point connections provide at least two 
linear constraints. All of the constraints of this type and all of 
the associated def ormational equations may be gathered together. The 
techniques of Section 3.6 may then be applied to any subset (including 
all) of those constraints and reiterated as required. As noted, when 
the coefficient matrices can he approximated as constant, the techni- 
ques outlined in Sections 2.3 and 2.6 may be applied to obtain modes 
of the coupled substructures. 

The procedures outlines in Section 3.3 - 3.5 may be applied to 
incorporate the nominal constraints for all of the point connections. 

As noted, it is generally possible to eliminate the constraint forces, 
at the expense of complicating the equations. However for rigid con- 
nections or for the variables along the axis in single axis connections, 
it is generally desirable to eliminate them. So, in general, some of 
the constraint forces will have been eliminated. In this event the 
number of equations is reduced, the number of variables is reduced, the 
constraint equations need not be retained, and the constraint forces 
do not appear. Those constraint forces which are not eliminated are 
treated in the next section (4.7). 

In the event the def ormational constraint motion is unknown and 
governed by a force-deflection relationship compatible with the form 
of the equations, the relative motion variable can be included in the 
deformation variables, and the governing relationship in the assembled 
equations. Analysis may then proceed as above. 
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4.7 Pfaffian Constraints 

Those single axis connections for which an election was made not 
to eliminate constraint forces, may be treated in the same way as 
Pfaffian constraints. The procedure is outlined in Section 3.2. The 
effect of this class of constraint is to augment the system equations 
by the constraint equations, and to augment the system variables by 
the constraint forces (Lagrange multipliers). As noted, unknown motions 
may be treated in a similar fashion. 

4.8 Equations for Structure 

After applying all these constraints, the number of frame variables 
is equal to the number of rigid body degrees of freedom. Similarly 
the number of rigid body modes remaining in the deformation equations 
is the number of rigid body degrees of freedom unless the option to 
reflect total relative motion in a frame variable (Sections 3.3 and 3.4) 
was exercised. For each of these cases where the nominal motion is 
arbitrary, but the relative deformation is constrained, the number of 
rigid body def ormational modes is reduced by one. The rigid body 
deformational modes can now be eliminated by selecting a mean motion 
or principal axis frame for the composite substructure, by fixing 
some degrees of freedom, or they can be retained. The method outlined 
in Section 2.5 can be used to integrate the coupled system in appro- 
priate cases, if all of the rigid body degrees of freedom are retained 
in the deformations. 

At this point, any structural elements which do not fit the 
procedure outlined above (e.g. non-Pfaffian constraints) can be incor- 
porated in the equations governing the structure. 
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4.9 System Synthesis 

The complete set of structural equations may now be combined with 
the equations governing the behavior of the other dynamical subsystems. 
To whatever extent the equations for the other subsystems share the 
form of the structural equations, the same techniques may apply. In 
particular, any subset of constant coefficient linear equations may be 
subjected to modal analysis. The complete set of system equations is 
now available for the intended analysis (e.g. numerical integration to 
determining system behavior). 
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APPENDIX A 


MATHEMATICS 

To preserve the continuity of the main body of this presentation, 
the mathematical definitions have been placed in this appendix. The 
material in this appendix is not new, and has not been originated by 
the author; rather it has been abstracted from the sources noted in 
the bibliography, particularly from Bishop and Goldberg's text "Tensor 
Analysis on Manifolds." The material is reproduced here primarily to 
establish the nomenclature, and secondarily, to review some of the more 
important results. For a more complete study, the noted sources are 
recommended. 

Throughout this appendix, the most fundamental aspects of math- 
ematics are applied implicitly. Some of these fundamental concepts 
and signs are noted below. Mathematics consists of mathematical 
objects and relations . The term set is synonymous with mathematical 
object. Mathematics is written in terms of complicated assemblies of 
fundamental signs and letters . A letter represents a totally indeter- 
minate mathematical object. The fundamental rules governing the use 
of the fundamental signs and letters are called axioms . The true 
relations which may be logically deduced from the axioms are called 
theorems . Additional signs are introduced to abbreviate complex 
assemblies. 

The fundamental logical signs may be taken as "or" and "not". If 
R and S are relations, the assemblies "R or S" and "not R" are rela- 
tions called the logical disjunction and negation . If at least one of 
R and S is true, then "R or S" is true. The negation of a true relation 
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is false . The derived logical signs: "and", and yield 

relations: "R and S", "R^^S" and "R^z^S", called logical conjunction , 
logical implication and logical equivalence ; and defined to be "not 
[(not R) or (not S)]", "S or (not R)" and "(R=>S) and (S=>R)", 
respectively. IE R is a relation, A is a mathematical object and x 
is a letter, then ,, (a|x)R m is the relation obtained by substituting 
A for x in R. If "(A|x)R" is true then A satisfies the relation R. 

The derived logical signs: "3" and "V" are called the existential 
and universal quantifiers and are read "there exists" and "for all". 
They yield relations: "(] x)R" and "(V x)R"; the first is defined by 
the axiom 

(A|x)R => (3 x)R. (A. 0-1) 

The second abbreviates "not [ ( 3 x) (not R)]". 

The fundamental mathematical signs may be taken as: "=" and " 6 " 

called the sign of equality and the sign of membership . The rules 
governing the sign of equality are 

( V x) (x = x) t 

( V x) ( V y) [ (x=y)0 (y=x) ] , 

( V x) ( V y) ( V z) [ (x=y and y=z) (x=z) ] and 

(u=v)=^[(ujx)R<^(v|x)R] . (A. 0-2) 

The rule governing the sign of membership is 

(A=B) <4^ [ (x € A)<^(x ( B)]. (A. 0-3) 

A.l Set Theory 

The fundamental concepts and symbols used in set theory are 
outlined in the introduction to this appendix. The symbol "C”, called 
the inclusion sign, yields a relation "A C B" which is an abbreviation 
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for "(Vx)[(x £ A)=^(x £ B)]". The relation is read "A is a subset 
of B" and is also written "B O A". The relation satisfies 
[(AC B) and (B C C)]=»(A C C) , and 

(A = B) [ (A C B) and (B C A)], (A. 1-1) 

If R{x} is a relation with x a letter (variable), then for every 
set, B, 3a unique A C B where 

x ( A<^(x £ B and R{x) is true). (A. 1-2) 

It is conventional to indicate this set by j x | R (x}j . The braces 
conventionally indicate a set and the vertical bar may be read "such 
that". In other uses, "such that" may be denoted " The complement 

of A in B (where A C B) or the set theoretic difference of B and A is 
denoted 

B - A . = jx|(x £ B) and (x i A)j. (A. 1-3) 

The set "0" = A - A is called the null or empty set. For any two sets, 
A and B, the intersection and union of A and B are denoted "A PI B" and 
"A U B" and are defined as 

A O B = x| (x £ A) and (x £ B) J and 

A U B = x | (x £ A) or (x £ B) J . (A. 1-4) 

N. N 

The notations U A, " 0 A ", " U A " and "Pi A " are used to 
i=l i=l 'aej a a€J a 

indicate unions and intersections of families of . sets. 

If x and y are mathematical objects, the mathematical object 

"(x, y)" is called an ordered pair . More generally, (x^,...,x ) is 

called an ordered n- tuple . The ordered pair satisfies 

[(w,x) = (y, z) [ (w=y) and (x=z)J. (A. 1-5) 

The cartesian product of two sets, A and B, is denoted "A x b" and 

defined as 


131 



(A. 1-6) 


A x B - j(x,y)|(x £ A) and (y £ B) 

More generally, x A^x.^xA^" is the set of ordered n-tuples. If 

two adjacent sets are the same, exponential notation may be used (e.g* 

2 

A x A x B = A xb). 

A set consisting only of ordered pairs (i.e. a subset of a 
cartesian product) is called a graph * A function , "f: A *> B", is a 
set. A, called the domain , a set, B, called the range and a graph, 

G,CA x g, satisfying V x £ A, 3 exactly one y £ B 9 (x,y) € G. The 
unique y associated with any x is denoted "f(x) n » The image , "f(A) n , 
is defined by 

£ (A) - jf(x)|x € A 1 . (A. 1-7) 

When f (A) = B, the function is called onto , otherwise it is into , If 
V y £ f (A) , 3 only one (x £ A) 9 (f(x) = y) , the function is one to one ; 
which is abbreviated "l-l". The function f f :C -*■ B, C C A, 9 f T (x) = 
f(x), Vx £ C is called a restriction of f and is denoted "f|„". If 

f: A B and g: C -* D, then ”g o f : E ^ D" is the function obtained by 

following f by g (i.e. g[f(x)]). The domain is E = fxj (x € A) and 
(f (x) £ C)j. The function M g o f" is called the composition of g and 
f. In the event E = 0, the function is the null function 0: 0 -+ D. 

If f: A -»■ B is 1-1 and onto, 3 a unique function called the inverse 

of f, denoted "f B -> A", ? f ^ o f (x) = x and f o f ^(y) = y. 

The algebraic development of the number system and its properties 
will not be pursued here. Notationally , the natural numbers will be 
denoted M N", the integer M I M , the rational "Z", the real "R ,r and the 
complex "C". 
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A. 2 Metric Spaces 

A, function,, p: A x A + R, is a metric if Vx,y,zfA, 

0 < p(x,y) < 

p(x,y) = 0 O x=y, 
p(x,y) = p(y,x) and 

p (x,y) < p(x,z) + p(z,y). (A. 2-1) 

Two metrics, p 1 , p 2 *. A x A + R, are strongly equivalent if 3 c 2 > 0, 
9 V x,y £ A, 

|^1 ( X, y> - c 2 p 2^ X,y -*] and [ p 2^ x,y ^ - C 1 p l^ x,y l] * (A. 2-2) 

Two metrics are equivalent if they give rise to the same "topology" 

(A. 3). A strongly equivalent metric is equivalent. A set with a metric 
is called a metric Space . The open ball centered at x with radius r, is 
defined by 

B(x,r) = {y|(y£A) and [p(x,y) < r]}. (A. 2-3) 


A bounded set is one which is contained in some open ball. 

A. 3 Topological Spaces 

A collection, t, of subsets of a set, A, is a topology in A if 


B.fcT, i=l, ...» N => O B.£T and 
1 i=l 1 

B^£ t , oi £ J => U B a £T (A. 3-1) 

where J may be finite, countable, or uncountable. Any set, B£t, is 
called an open set , A set with a topology is called a topological 
space. A topology constructed from a metric is called a metric topology . 
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If A and B are topological spaces, the topology whose open sets are 

cartesian products of open sets in A and B is called the product 

topology . The formal statement "A with T is a topological space" is 

customarily abbreviated as "A is a topology." 

A neighborhood of a point is defined as an open set containing the 

point. A closed set is the complement of an open set. For any subset, 

9, the union of all open sets, C^jCB, is called the interior of B and 

is denoted "a 0 "; the intersection of all closed sets, C , D B is called 

the closure of B and is denoted "B the set B - B° is called the 

boundary of B and is denoted "3B". A topological space, A, is connected 

if 0 and A are the only sets which are both open and closed. An open 

connected set is called a region . A function, f: A ■+ B, from a 

topology. A, into a topology, B, is continuous if for every open 
-1 

C C B, f (C) is open in A. 

Two topologies, A and B, are homeomorphic , and a function, f: A -+• B, 

-1 

is a homeomorphism, if f is 1-1, onto and both f; A B and f : B -*• A 

are continuous. This is the natural notion of topological equivalence. 

A. 4 Manifolds 

A chart at a point, p, 6 A, a topological space, is a homeomorphism, 
p: U *+ R^, mapping U, a neighborhood of p, onto an open subset of R^. 

The dimension of the chart is d. The standard Cartesian chart on R is 

d. d 1 2 

the identity map, R -* R , which is defined by (u^Cp), u (p), . .., 

d, ... .12 dx , i. 1 2 dv i 

u (p)), where p«* (a , a , a) and u (a , a , . . . , a ) =* a . The 

function u*: R is called the i-th projection . A chart is also 

called a coordinate map . The entries in the chart, p, are defined by 
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u o pand are called the coordinate functions , usually denoted by 

tt " i 2 d " 

x . Collectively the coordinate functions (x , x , ...,x ) are 

called the coordinates or coordinate system at the point. 

A function, f : V + R, V open in R^, is continuous to order 

" oo" 

infinity, denoted C , if f has continuous partial derivatives of all 

0 OO 

orders and types. A function, <p ; V R , is a C map if its entries, 

i oo d e 

uo<p, i - 1 , e are C . Two charts, p : U R , and T : V + R , 

OO 

on a topological space, are C —related if d = e and either U HV = 0 
—1 —1 °° 

or both p o t and Top are C maps. A collection of charts, 

{p : U ■+ R la 61, U CA), such that {U la€l} is a covering of the 
a a a a| - 

OO 

topological space. A, (i.e., AC U U ) is called an atlas . A C 

a € I 

OO 

atlas Is one for which every pair of charts is C - related. A chart is 

OO OO 

admissible to a C atlas if it is C - related to every chart in the 
atlas . A topological manifold is a separable Hausdorff space 93 a 
d-dimensional chart at every point. The dimension of the manifold is 
d. (A metric topology is a Hausdorff space, and the metric topologies 

OO 

in this presentation are separable.) A C manifold is a topological 


manifold with all of the admissible charts to some C atlas. 

If M and N are manifolds of dimension d and e, then the product 
manifold is the manifold of dimension d+e obtained from the product 
topology with an atlas obtained from the products of charts from the 
atlases of M and N. The product of charts, p: U -*■ R^ and t : 

V + R , is (p,x): U x v ■> R , where 

[Cy,T>] (m,n) - [p(m), t(n)j. (A. 4-1) 
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oo 

The definition of a C map may be extended to' include a function, 

oo J 

F : M ->■ N, where M and N are C manifolds. Let y^j u -+ R and 

£ OO 

M 2 : V + R be C charts on M and N. The y^ - \i^ coordinate expression 

for F is 

y 2 ° F ° : o F _1 (V) -> R e . (A. 4-2) 

If all such expressions, for all admissible charts, y^ and are C°° 

OO 

Cartesian maps, the map F is called a C map . (It is sufficient that 
the coordinate expressions for one atlas in M and one in N are C°°.) A 

OO 

diffeomorphism is a 1-1, onto, C map, F : M -* N, such that the inverse 
—1 00 

map, F : N -+■ M, is also C ; M and N are called dif feomorphic . This is 

the natural notion of equivalence between manifolds. 

0° 

If there is a 1-1, C map, F : M -+■ N, 9 at every m£M] a neighbor- 
hood, U of m, and 3 a chart of N at F(m), y : V -> R e , y = (y 1 ,..., y e ) 

9 x » y* o F ly, i = 1, . . . , d are coordinates on U for M, the map, F, is 

called an imbedding. A submanifold of N is a subset, F(M), where F: 

M •* N is an imbedding, provided with the manifold structure for which 
F: M F(M) is a diffeomorphism. An open submanifold is an open sub- 

set of a manifold, with the manifold structure obtained by restricting 
the topology and the coordinate maps to the subset. 

A. 5 Tangents 

A differentiable curve , y : [a,b] -* M, is a map of an interval of 

real numbers, [a,b], into a manifold, M, ? 3 an extension to an open 

interval, y : (a-c, b+c) -*• M (c > 0 and y - y on [a, b]), which is a 

OO 

C map. 

OO 

The set of all real valued C functions, F ; U R, where U is a 

l» oo 11 

neighborhood of a point, m, € M, is denoted F (m) . A tangent at m is 
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a derivation of F (m) , which is a function, t : F (m) R , ) 

Vf,g€F (m) and Va,b£R 

t(af + bg) = at(f) + bt(g) and 

t(fg) * t(f) g (m) + f(m) t(g). (A. 5-1) 

03 

If Y is a C curve in H ? y{c) = m, then the tangent to Y at c » 

' » 00 
Y*( c ) » is defined by, Vf £ F (m). 


[Tf*(c)I(f) - (O. 


(A. 5-2) 


The set of all tangents at a point, m, is called the tangent space 

tl II 

at m and is denoted M 

ra 

OO 

A real valued function at m, f, € F (m) , has a coordinate expres- 

1 d “1 d 

sion, f — g(x, x ), where g - f o ]J : U R, U open in R . The 

partial derivate of f with respect to a coordinate, x'*', is defined as 

(A. 5-3) 


3 3 , N 3(f o \i l ) 

“ T (f) = r (g) 0 V = 7^ o \1 . 


3x 


Su 1 


3u‘ 


i" 


The operators 3/3x are called the coordinate vector fields . If 

3 

application of — ~ is followed by evaluation at m, the result is a 
3x n ^ ii 

tangent at m, which is denoted — r-(m) , defined by 

3x 


(m)J(f) = (f)J(m). (A. 5-4) 

Specifically, it is the tangent to y^ at m, where is the i-th 
coordinate curve through ra defined by 

Y, (u) = p 1 [x 1 (m) , . . . , x 1 1 (m), u,x 1+1 (m) , . . . , x d (m)]. 

(A. 5-5) 

These tangents form a ’’basis” (A. 6) so that for any a 1 £ R, i«l, d 
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t * £ a* — (m) 
1=1 9x 


(A. 5-6) 


is a tangent. Conversely, for any tangent, t,fM , 3 unique a 9 equa- 

m 

i ' 

tion (A. 5-6) holds. These a are called the components of t with 

i 

respect to the coordinates x « 

It is convenient to introduce an abbreviation for sums called the 
Einstein summation convention . The appearance of an index exactly twice 
in a product (usually once up and once down) implies summation over the 
range of the index. For example, if a 1 and b- are defined for i=l,... f 
N, then M 


aV 


£ i, 


1=1 a \- 


(A. 5-7) 


Thus the expression for any tangent may be written as 
i 9 


t = a 


T (m) . 


(A. 5-8) 


9x 


The a may be determined by noting that 

Ux* J } 

1 


(x^) = a 1 PL (m)j(x^ 


= a 


- a J 


6 , 


(A. 5-9) 

The symbol "d'j" is called the Kronecker delta and is defined to be one 


If i= j and zero if i ^ j . The contraction is 


' S l + 6 2 + S 3 " 3 - 


(A. 5-10) 


The tangents are algebraically "vectors" (A. 6) and the set of 
tangents at a point forms a "vector space" (A. 6). In addition, the 
components with respect to the coordinates satisfy the usual tensor 
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then the 


t r an s f o rma t ion relations. If (x^) and (y d ) are coordinates 

g i 

tangent — t (m) has an expression in terms of the x , 

9y 

““T (m) = (m) , 

9y 9x^ 


(A. 5-11) 


where 


So 


j . { jx 9x^ , . 
r = t(x J ) = — - (ra) . 

9y 


9 / \ Sx- 3 , v 9 { . 

— r (m) = — r W — r (m) . 

9y 9y 9x^ 


(A. 5-12) 


(A. 5-13) 


Thus for any tangent, 
i 9 


t = a — r- (m) = b 1 — — (m) 
3x 9y 


1 _ J ^2L 

3y 


a = b 


j 


(m) , 


(A. 5-14) 


which is the customary tensor transform. 

The set of all tangent spaces at all points of a manifold, 

M, is called the tangent bundle , and is denoted "TM". If p: M — ► N is 

oo 

a C map, then the tangents in M are mapped to the tangents of N by 

P*: TM'** TN, the differential of p . In terms of coordinates, x^,i=l,..., 
ot 

d, at m and y , a = 1 , e at n = p(m), if p has a coordinate expres- 

sion, y°op = f a (x 1 , ...» x d ) , then a tangent, t,Oi , t = a i — ~ (m) is 

. m 9x L 

mapped to the tangent p*(t)€N , (n=p(m)), where 


n 


y*(t) = b a (n) and b a = 
9y a 


= a 1 (y a °V>)J ( m > 


(A. 5-15) 


9 a 

The array — 7 (y op) is called the Jacobian of p with respect to the 
9X 1 

i ot 

coordinates {x } and {y }. In the special case of a real valued func- 
tion, f: M ->■ R, 


f*(t) = t(f) (c), where c = f(m). 


(A. 5-16) 
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The component, t(f), is redefined as the differential of f and the 
notation is changed to 


[df ] (t) = t(f). 


If lx 1 are the coordinates on M, 

df = — — r (f) dx 1 

ax 1 


and 


w -‘ 1 ( i >) = 


(A. 5-17) 


(A. 5-18) 


(A. 5-19) 


A vector field, X, on E C M is a function assigning to each mf E 

QO 

a vector, X(m), € M . It is C if the components of X with respect to 

i 00 

every coordinate system, x , are C functions. A curve, y, is an inte- 
gral curve of a vector field, X, defined on ECM if the range of y is 
in E and for every s in the domain of y 

Y*(s) = X[y(b)J. (A. 5-20) 


A. 6 Vector Spaces 

A commutative group is a set, X, with a law of composition, 

(x,y) _► x+y, Vx,yCX satisfying, Vx,y,z t X, 
x + (y + z) = (x + y) + z, 

X + y = y + x, 

30?x + 0 = 0 + x = x and 

3 (-x) ? x + (~x) = (-x) + x = 0. (A. 6-1) 

A r ing is a commutating group with a law of composition, (x,y) xy, 
Vx,y€X, satisfying, Vx,y,z£X, 
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x(yz) = (xy)z, 

3 1 9 lx = xl = x, 

x(y + z) = xy + xz and (x + y)z = xz + yz. (A. 6-2) 

A vector space is a set, V, whose elements are called vectors , 
a ring, S, whose elements are called scalars , and two operations, 
vector addition , (x,y) ->•' x + y, V x, y ( V and scalar multiplication , 
(A,x) ■> Ax, Vx (V, AfcS such that V with addition is a commutative 
group and Vx,y (V and A,p € s 
A(yx) = (Ap)x, 
lx = x, 

(A + p)x = Ax + ]ix and 

A(x + y) = Ax + Ay. (A. 6-3) 


The ring of scalars will be taken to be R, the real numbers, unless 
explicitly noted otherwise. 

i i 

A sum of products of scalars, a , and vectors, v^, i.e., a v^ is 
called a linear combination . If all of the a 1 = 0 it is trivial. 
otherwise it is nontrivial . A finite set of vectors {v^j is linearly 
dependent 3 a nontrivial, null linear combination, a 1 ^ = 0; other- 
wise the set is linearly independent . A non-empty subset, U, of vector 
space is a subspace if it is c losed under addition and scalar multi- 
plication (e.g., V u,vf U,a6S, u + v£U and au£U). The minimal sub- 
space containing a subset, T, is called the subspace spanned by T. 

A linearly independent set of vectors, T, spanning a space, V, is 
called a basis for V. It can be shown that every basis for a space 
has the same number of elements, d, the dimension of the space. 
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Any vector, v£V may be represented as a unique linear combination 
of basis vectors for any basis (v i >, v - a v^. If (e^ and {f^} are 
two bases of V then each basis vector, e^, has an expression in terms 
of the {f } , 

e. = a^ f , . - (A. 6-4) 

i i J 

Similarly, 

f. « e.. (A. 6-5) 

i i j 

If the scalars, a^, are arranged in an array with j constant on rows 
and i constant on columns, the array is called a matrix and is denoted 

n - M 

(a^) . Specifically, it is the matrix of change of basis from 

{e . ) to if.}- Since 
1 J 

e = a*! (b k e,) , (A. 6-6) 

i l j k r 

it may be concluded that 

a.] = 6^. (A. 6-7) 

1 j i 

Similarly, 

b] a k = 6 k . (A. 6-8) 

1 J i 

Mat rices, (a^) and (bt^), satisfying these relationships are called 
inverses . 

A. 7 Tensor Spaces 

If V and W are vector spaces, f: V W and Vv^, v 2 £ V and a £ R, 

f(v^ + v 2 ) - f(v^) + f(v 2 ) and 


f (a v 1 ) = a f (v^) 


(A. 7-1) 



then f is a linear function . The vector spaces V and W are isomorphic 

and the function f is an isomorphism if f is a 1-1, onto linear func- 

_1 

tlon. The set of points mapping into the zero of W, f (0), is called 
the null space of f. The set of all linear f: V W, {f |f: V^W and 
f linear}, is called the space of linear functions from V into W and is 
denoted ”L(V,W>”. Defining [f + g](v) = f(v) + g (v) and [af](v) = a f 
(v) as vector addition and scalar multiplication, L(V,W) is a vector 
space. 


Just as the selection of a basis provides a coordinatization of a 
vector space, selection of bases, ie^}, i = 1, * . • , d^ on V and {e^}, 
a = 1, ..., 6-2 on W (where d^ and are the dimensions of V and W) 
provides a coordiatization of linear functions as x ^2 matr i ces > 


( f “)’ 


Hep = f“ (?„>. 


A basis for L(V,W) is {E^}, where 
e s= <5^ e 

L e i i 


(A. 7-2) 


(A. 7-3) 


For any f, 


£ " E 6 


The matrix of E^ is 6^ 


(A. 7-4) 


The set of scalar-valued linear functions, L(V,R), is called the 

dual space of V and is denoted n V* " For any given basis of V, {e^}, 

there is a unique basis of V , such that is 

called the dual basis to {e,}. The dual basis to the coordinate 
! 

, is a basis for the space of differentials . 


vector field. 


Sx 1 


r (mV 
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idx^On)} , (i.e. dx^ 3. = 6^). If if.) is another basis of V and 

J J * 

its dual, where f. = aj e, # and e. = f. and (b“! J is the inverse of 
iij* lij \ i / 

j^a^J , then <p and = b^ . The dual space to V (i.e. 

^ £ 

(V ) ) is naturally ismorphic to V. The two will not be distinguished 

A A 

and (V ) will be denoted V. 

If V^, i=l, ...» r, and W are vector spaces, f : x ... x -»• W, 

and Vv., v. 6 V.. and V a, a £ R, 

iii ’ 

av i + a v r ) = a f(v v ±9 . . . , v r ) 

+ a f(v 1? v ± , . v r ), 

(A. 7-5) 

then f is a multi-linear , specifically an r-linear function. The 

II 

space of r-linear functions is denoted "L(V^,..., V^; W) . The scalar 

* 

valued multi-linear functions with variables all in V or V are called 

tensors over V and the vector spaces they form are called tensor spaces 

* 

over V. The number of variables taken from V is called the contra - 

variant degree, and from V the covariant degree . The functions, f : 

* * * 

V x tt ,xV x Vx ... xV -*■ R (V r times, V s times) form a tensor space 

denoted M T r (V)", 
s 

T^(V) = V©...©V© V*® ...©V* (A. 7-6) 


/v 

(V r times, V s times). They are called tensors of type (r,s) and have 

contravariant degree r and covariant degree s. 

A tensor of type . (0,0) is called a scalar , Tq = R; of type (1, 0) 

1 * 

is called a contravariant vector, T^ = V = L(V ,R); of type (0,1) is 

A ^ _ 

called a covariant vector , - y = L(V,R)» The symbol, ” (x) n , 

IT t 

denotes the tensor product which assigns to tensors A ( T and B £ T 
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a tensor A <*) B € The tensor, A © B ; (V*) 1 L x (V) ST%1 -*• R, 

is defined by 

A©B (x 1 ,..., T r+t , v g+u ) 


r+t 


s+ii 


= ACT 1 , . .., T r , V g ) B(T r+1 , ..., T r+t , v 


v ) 
s+l’ S+u' 


(A. 7-7) 


Selection of a basis, {e^}, and its dual, {e^}, provides a coordina- 


tization of a tensor, A,6T g , as a mult i=dimenslonal array; 


i i 


A = 


A. 1 , ... r e i 0 ...© e i ©e^l @.,, 0 e^s 
J 1 -*s 1 1 r 


where 


in » * « r 


'V*d. ( ’ 


, £ r, e 


, • • • > e . 

1 J s/ 


(A. 7-8) 


(A. 7-9) 


If the basis is changed to {f^} with an associated dual {y 1 }, where 

f. = a") e, and v* - b 1 £^ , then the coordinate expression for A in the 
i x X j 

{f .} basis is 
i 


VV" m r - b"l ...b m r a J l ... ah A 1 !"^ 
f i, in, n e j - . . • i 

n....n 1 r 1 s J 1 

I s 


(A. 7-10) 


A tensor is called symmetric if it is both covariant and contra- 
variant symmetric, a tensor is (co- or contra-) variant symmetric if 
it is symmetric in every pair of (co- or contra-) variant indices. It 
is symmetric in the p-th and q-th (co- or contra-) variant indices if 
the components with respect to every basis are unchanged when the 
indices are interchanged. If the interchanging changes only the sign, 
the tensor is skew-symmetric in the p-th and q-th (co- or contra-) 
variant indices. It can be shown that symmetry with respect to one 
basis implies symmetry with respect to every basis. 
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A tensor valued function of tensors is frequently discribed in 
terras of the components with respect to a certain basis. If the result 
is independent of the basis, the function is called an invariant . A 
contraction is an example of an invariant. A contraction assigns to a 
tensor of type (r,s) another tensor of type (r-1, s-1) . Formally, a 
contraction in the p-th contravariant and q-th covariant index is 
defined by 


n i. * * . x -* . i . . . * 1 1 k 1 

B.l . r-l = A.l , p-1 p 

J r*° s -i J r*°q-l k \ 


. i 
•j 


r-1 
s-1 • 


(A. 7-11) 


The scalar product of a vector and a dual vector is an example of a 

i i 

contraction. If v =* v e. and t - T. C then 

i J 

V©T - v T. e. © (A, 7-12) 

•> J 

r 

and the contraction is v x. , There are several alternative inter- 

k 

* 

pretations. Recalling V = L(V,R) the scalar product may be thought of 

•A? 

as a function, T (V) , or with V = L(V , R) as a function V of x, or as a 

* 

bilinear function < , > ; V x V -► R defined by 

<v , T> = T(v). (A. 7-13) 

An invariant may be linear or multi-linear. A linear invariant of the 
p-fold tensor product of a variable with itself is called an invariant 
of degree p . 

A quadratic form on V is an invariant of degree 2 with variable 
in V . A tensor of type (0,2) is called a bilinear form , b : V * V -► R, 
and has a coordinate expression 

b = b^ e i ®e^. (A. 7-14) 
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For every quadratic form, q(v), there is an associated symmetric bi- 
linear form, b, defined by 

b(v, w) = j j^q(v-bw) - q(v) - q(w)^j . (A. 7-15) 

Conversely for every symmetric bilinear form, b, there is an associated 
quadratic form, q, defined by 

q(v) = b(v, v). (A. 7-16) 

A bilinear form is non-degenerate if for every v€V, v ^ 0, 3 w(V?b 

(v,w) ± 0, or equivalently, if the matrix (b,^) is non-singular. A 

quadratic form is positive or negative semidefinite if V^v ^ 0, q(v)> 0 

or q(v) < 0, respectively; if q(v) ^ 0 for v ^ 0 the form is (positive 

or negative) definite . A definite form is non-degenerate. 

A non-degenerate symmetric form (usually positive definite) is 

called an inner product . A vector space with a positive definite inner 

product is called a Hilbert space and implies a topological structure. 

Two vectors, v,w, are orthogonal with respect to b if b(v,w) - 0. If 

b (v,v) = 0 then v is called a null vector of b . A basis, {e^}, of V, 

is orthonormal with respect to b if for i ^ j , e_. ) = 0 and 

b(e e ) . no sum, is one of the values + 1, 0, - 1. The bases, 9 (m) , 
a* a 

in general, may not even be orthogonal. However, an orthonormal basis 
always exists, if the form is definite, the Gram— Schmidt process 
yields an orthonormal basis. If, as in the usual orthogonal curvi- 
linear coordinates, the basis is orthogonal (and definite) but not 

-1 1 

normal, the change of basis, = h j^On), yields an orthonormal 
basis. The h^ are the Lame coefficients defined by 
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3x^ 3x 


(A. 7-17) 





, ,-l k, j 
and h , h j 

j 1 



A metric on a manifold is a bilinear form, b, ’which is defined and 
non-degenerate at every point of the manifold. For the usual coordin- 
ates on E, S or T, 

b = 6., dx"*" (*) dx 1 ^ (A. 7-18) 

ij 

(i, j-1 , 2, 3, 4 for E, i, j - 1, 2, 3 for S, and i, j = 1 for T) . 

This metric is Reimanian (positive definite) and the pair (M, b) is 
called a Reimanian manifold. The inner product or dot product is 
frequently denoted by "< , >" or n • " ; 

b(v,w) = <v, w> - v • w. (A. 7-19) 

The use of the dot product is extended to tensors (e.g. A • B) by 
interpreting it as acting on the last contravariant index of the pre- 
ceding tensor, and the first of the following tensor. There are three 
interpretations of the action of the fundamental bilinear form repres- 
ented by the dot: 1) lowering the last contravariant index of the pre- 

ceding tensor and contracting, 2) lowering the first contravariant 
index of the following tensor and contracting, or 3) directly assign- 
ing a tensor of contravariant degree two less. As an example, 


N * 



a$ d 
3x 


a 


© 


3x 


= n y s 

oe 


aP 0-\ © — g ®dx 6 ® dx £ . 

8x y 3x U 3x P 


(A. 7-20) 


Using the first interpretation, 
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(A. 7-21) 


N * S = n s 
— = e 


a$ 3 


3x 


a 


© 


(x) dx* 


3x h 


- n s 

a 


a£ 3 


3x 


3 ’ 


The length of a tangent on a Riemanian manifold is 


v] 1 - <v,v>" . 


(A. 7-22) 


The angle between vectors is the value of 8, 0 < 0 < tt, such that 


/ns < v, w> 
cos (0) = 1 


(A. 7-23) 


I M I I |w| | 

The length of a curve , y : [a f b] •> M, denoted "|y|", is the integral 

of the lengths of the tangent vectors, 

.b 

i I S/ f j-\ I I J*. 

(A. 7-24) 


Y = 


-f I ly*c 


t) dt. 


The distance between points is the greatest lower bound of the 
lengths of all curves. The Riemanian metric, b, is consistent with 
the previous definitions of distance and duration, and straight lines 
are minimum distance curves. 

The usual coordinates make E,S and T essentially equivalent to 
4 3 d 

R , R and R. If R had been taken as the models of the physical 
system, or if the metric is carried over to R^, then b would be the 
standard flat metric and (R^,b) would be ordinary Euclidean d-space . 

A. 8 Tensor Fields 

For each type, (r,s), of tensor, and each m£M,M. a manifold, there 
is a corresponding tensor space M ^ ms " over M^, the tangent space. 

For any (r,s) the union of these tensor spaces Vm£M is called the 
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bundle of tensors of type (r,s) over M, and is denoted "T M," In par- 
ticular the tangent bundle is TM = T^M, and the bundle of differentials 
is ' T^M. A tensor field, T, of type (r,s) is a function, T : E ->■ T*M 
where E C M. In particular, X: E->TM is a vector field, and the dif- 

do 0 

ferential of a C function, f, df : E**T^M is a tensor field of type 

( 0 , 1 ). 

oo r I"t 

A vector field is "C " for every coordinate system lx ) , the 

i i i 00 

components of X with respect to the x , X = X (x ) , are C functions. 
The components of T £ T*M, with respect to the coordinates {x^}, are 
the d real valued functions. 

T 1 !-" 1 ! = T^dx^l, . , dx i r , — - 4 — ). (A. 8-1) 

j . . . j \ dx^l 9x J s f 

J 1 J s 

00 DO 

The tensor field, T, is C " if its components are C functions. A 

oo 

tensor field of type (0,1) that is C is called a 1-form (or Pfaf f ian 
form ) . 

OO 

If ip : M + N is a C map, then vector fields, X and Y, on M and 
N are co related <=> V m in the domain of X 

tf*[X(m)] = Y [ (/> (m) ] . (A. 8-2) 

The definition is extended to tensor products by 1 

i P* m ( A ® B ) = (A)] ® [ & * m (B) ] • (A. 8-3) 

Then two contravariant tensor fields, S, on M, and T, on N, are 
(p -related if Vm£ domain of S. 

[S (m) ] - T [ <p (m) ] . (A. 8-4) 

For covariant tensors, 

i/?* m (T[ (/>(m)]) = S(m), (A. 8-5) 
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(A. 8-6) 


it ■ft it 

where </> : N ip{ m) -+ M m is defined by 


m 


< 'P* m (v),T> = <v, fm (T)>. 


In terms of coordinates if 


y o tp = F (x , . . . ,x ) 


and 

T * yG) d y a @ dy ^» 

Then 

<A*> = T (F 1 (x) , . . . , F 6 (x)) —^rdx 1 © dx j 

Sx 1 8x J 


(A. 8-7) 


A, 9 The Hodge Star Operator 

The space of skew- symmetric tensors of type (r,o) is denoted 
by a wedge, "A V," where V is the associated vector space. The dimen- 
sion of A r V is the binomial coefficient M ^ d where d is the dimen- 
sion of V. If Is a permutation of i^,..., i^, it may be 

obtained in any of a number of ways by transposing pairs of indices. 

The number of transpositions required for a permutation is odd or even 
and the corresponding sign of the permutation is -1 or 1, resectively. 
For a permutation, it, this is denoted "SGN(tt)." The component of a 

■? n J* # 

tensor, A J l'*' J r* is the same as, or the negative of Al*V* 1 r if 

SGN(j 1 ...j r ) is +1 or -1, respectively. The alternating operator , 
r r 

A: T -► A V, assigns to each tensor, B, its skew-symmetric part, 

denoted ”B ". For 8-...G £V B is defined as 
a 1 r a 


B (0 
a 


1’ 



*'<F *?"(*■ %)■ 

(A. 9-1) 


where the sum runs over all r! permutations of (l,...,r). The 
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exterior product , denoted by ” A", is defined by 


A A B = (A © B) , (A. 9-2) 

a 

where A and B are skew- symmetric . The treatment of covariant tensors 
is analogous. 

An orientation of V is given by a non-zero element, 0, of 

d 3 

A V. For the usual coordinates on R 


0 = _L. A 8 A JL_ (A. 9-3) 

dx dx dx 

is taken. An ordered bais, (e^,..., e^), is in the orientation given 
by 0 if 


e^ A,,. A - a 9, a > 0, 

Since the binomial coefficients, ^ j and | j 


(A. 9-4) 

are the same, the 


it d —17 

spaces of skew-symmetric tensors, A V and A V, have the same 

dimension. The Hodge star operator is an isomorphism between the 

two spaces. Let e^) be an ordered orthonormal basis in the 

orientation. A typical basis element of A V is a ...A e^ . Let 

1 r 

jp..*, be chosen such that (i^,..., i^, r ) is an even 

permutation of (1,..., d). Then * is the linear transformation such 


that 


^ (e^ a • • . a ) ° 6, a ... a e, . (A. 9—5) 

1 r ~*1 

It may be noted that for an odd dimensional space (e.g., 3),* o * is the 
identity. The dot product may also be expressed in terms of the * 
operator, 

a • b - * (*a Ab). (A. 9-6) 
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If a and b are vectors of dimension three, the vector or cross 
product is defined to be 

a x b = *(a a b) . (A, 9-7) 

By definition, if a = a 01 - ■ — and b = b^ ^ 


9x 


a 


9x 


a* 


a A b = (a (x) b) 

( a , 
a b 


a 

3 3 


3 


a ® g 
3x a 9x P / a 


■) 


(a“ b B 

- a B b a ) -2- 

® - 

V 

f 9x a 

9 

2,3 3,2 

>\ 9 A 9 , / 

3,1 

. b -a b 

) 2 3 + \ 

9x 9x J v 

a b - 


9x 9x" 


dx 9x‘ 
(A. 9-8) 


So in terms of components, a x b is 


/ 2.3 3.2 \ 9 ./ 3,1 1.3 \ 9 . / 1,2 2,l\ 9 

* (a A b) = ( a b -a b 1 — ■£ + ( a b -a b 1 — j + I a b -a b j — 
■ / 9x v / 9x V / 9x 

(A. 9-9) 

This may be abbreviated by use of the Levi-Civita epsilon , e:^ , 
defined as: + 1 if (a,3,Y) is an even permutation of (1,2,3); -1 if 
it is an odd permutation; and zero otherwise. 


3 * 


a x b = a^b"^ — , 


(A. 9-10) 


The Levi-Civita epsilon may be used to express the determinant of 
a matrix of order three. 


DET ^ A ij \ e ijk a li a 2j a 3k 


e ijk a il a j 2 a k3* 


(A. 9-11) 
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The rule for row or column Interchange may be expressed as 


a, a J a J 

lp lq ir 


a . a , a . 

jp jq j r 


^p ^q a kr 


E ijk £ pqr |A 


(A. 9-12) 


The contraction on the first indices may be obtained by setting 

a m = ^mn* 



5 

<$ 




6 m 


6 _ „ 

6 

6^ 

°11 

lq 

lr 


22 

2q 

2r 


33 

3q 

3r 

< 5 ., 

6 . 

6, 

+ 


6. 

6. 

+ 

5.0 

6 ., 

6 

jl 

iq 

jr 


j2 

jq 

jr 


j3 

jq 

jr 


6 

<5 


6 

6 , 

6 


6 „ 

6 

6, 

°kl 

kq 

kr 


k2 

kq 

kr 


k3 

kq 

kr 


£ ijk £ iqr 


(A. 9-13) 

For non-trivial cases, q ^ r and the first column, will be identical 
to the second or third in two of the arrays. The third reduces to 


£ ijk E iqr = 


6 6 . 

jq i r 

6 . K 

kq kr 

6 . 6 - 6 . 6 . 

- jq kr jr kq 

A further contraction gives 

ijk ijr jj kr jr kj 

- 36, - 6, 

kr kr 


= 26 


(A. 9-14) 


kr . 

A final contraction gives 
E ijk £ ijk " 26 kk 


(A. 9-15) 


= 6 . 


(A. 9-16) 
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The results In Equations (A. 9— 11, — 12) yield some useful 
identities for direction cosine matrices (whose determinants are one), 
namely 

E ijk E pqr c ip c jq c kr 

= E pqr C pi C qj c kr’ (A.9-17) 

A. 10 Integration 

^ f f erential p-f orm is a C , covariant, skew-symmetric tensor 
field of degree p (type (0,p)). For the case of R 3 , the bases for 
forms are: 

0 - forms: {1} (i.e. the real valued functions) 

1 - forms: {dx, dy, dz] (i.e. the space of differentials (A. 7)) 

2 - forms: {dy a dz, dz a dx, dx a dy} 

3 - forms: {dx a dy a dz}. 

The Hodge star operator (A* 9) may be used as it was for contravariant 
tensor fields. 


The exterior derivative of a p-form, 0, is the p+1 form 
(i- 


d0 = d0 , . . x A dx^ 1 ... a dx lp , 

L i • m m X / 

1 P 


(A. 10-1) 


summed over all increasing sets of indices. 

3 

As an example, on R , if 
0 = fdx + gdy + hdz. 


then 
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d0 - df A dx + dg A dy + dh A dz 



(A. 10-2) 

Integration on a manifold is defined only for forms. Forms are 

oo co 

integrated over sums of oriented C p-cubes . AC p-cube is a map- 

P 

ping, a : U -+ M from a rectilinear p-cube in R , U, into the manifold, 

M. A rectilinear p-cube in R p , D, is, the set 

U « {(u 1 ,..., u P )| h 1 < u 1 < b 1 + c 1 , i=l, . . . , p}. (A. 10-3) 

The integration is performed by pulling back the differential form 
to R p (i.e., ot* (0 ) ) • A scalar valued function on R p is obtained by 
defining an inner product on p-forms such that the orientation, to, 
is unitary, 

<w, w> p = 1, (A. 10-4) 

and taking the inner product with the orientation, 

<a*(0), co> p . (A. 10-5) 

OO 

The integral of a p-form, 0, on an oriented (co) C p-cube, a: U •* M, 

is denoted 



(A. 10-6) 
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and is defined to be the associated Riemannian integral 


/ e -/ 


U < ot*(0),o) > p dy p . 


(A. 10-7) 


The generalization of Stokes 1 theorem to forms is as follows: 
let 0 be a (p-l)-form defined on a chain of p-cubes, C,(p>0), then 


A de = Ac 8 • 


(A. 10-8) 


A. 11 Variations 

A real valued function of functions is called a functional . 

(e.g., y < v, v > dy^ is a real valued function of the functions 
M 

v.) The variation of a functional, , is denoted 

"6J [h^ , . . . ,h^] n and is defined as the principal linear part of the 

increment . 

aj [h^,.,.,h^j = J[y-^ + y N + - y ^ * ... * 

: (a. li-i) 

is the increment of the functional, J, corresponding to the incre- 
ments 

h i (x 1 ,.., x d ) of the functions y (x 1 , . . . ,x d ) . 

The increments are arbitrary members of the noraed function space 
appropriate to the associated function y. (e, 


: ' 8 " A F (ipf) 


etc., suggests the space of 3-times continuously differentiable func- 
tions over R, D^(R).) If the increment can be expressed as 


157 



AJ[h 1 ,...,h n ] = . ..h^] + c( | | h j | ) 


(A. 11-2) 


e -+• 0 as | | h J | + 0; and (J) is a linear functional, then <p is called the 
principal linear part and is thus, by definition, the variation of the 
functional. A functional, (j), is linear if is continuous and 


Hh. 


ah^+ 


$h 


i» 


V 


= a(J>[h n 


» h . 


V 


+ 34>£h 1 , . . - >h ± , . . . ,1^] 


(A. 11-3) 


A natural norm for h is 

llhll - E l|h ± 

i=l 1 


(A. 11-4) 


and for an h. £ D ,(R) is 
i M 


MAX 


MAX 


|h t ! + ... + E ••• E _ 


R 


J l =1 V 1 R 


A i 

Sx . . . . 3x . 


|. (A. 11-5) 


*M 


Alternatively maximums could be used in place of sums. The term 
G ( | | h | | ) , if £ +.0 as | | h| | +0, can be written using little - 0 
notation as "o(||h||)," meaning 


L1M 1.1 h J 1 ) » o , 

I I h | | - 0 ||h|| 

or with big - 0 notation as "0( | |h| | 2 ) , meaning 3 m ? 
e||h|j < M||h|| 2 . 


(A. 11-6) 


(A. 11-7) 
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A. 12 Nab la 


The symbol "V u , called nabla or del , is used to represent the 
operation of the exterior derivative on differential forms (A. 10). 
The gradient is denoted n Vf n and is defined as df. The divergence 
is denoted "V_ » _V" and is defined as * d * V, where is the Hodge 
star operator (A. 9). The curl is denoted M _V * V n and is defined as 
* dV. 

In terms of the usual coordinates on three dimensional space, 
these have the usual forms 


,/r\ 3 f J . df J df , 

d(f) = — dx + dy + — d Z , 


(A. 12-1) 


_V . v = * d * (V^dx^ + + 


■ 2 3 3 1 1 2 

= * d (V^dx a dx + ^2 dx a dx + dx a dx ) 

^ v i 1 2 3 ^ V 2 2 3 i 

= * ( — _ d x a dx a dx + — t dx a dx a dx 


dx"* 


dv. 


3 3 1 

+ — - dx a dx a dx 

3x J 


dx 

> 

2 


= * 


fa 

lax 1 


dV r. 3V~ 

4- H 

2 3 

dx dx 


j dx 1 


2 3 

a dx a dx 


1 . 2 + 
dx dx dx 


a JH 


and 


(A. 12-2) 


V x V = 


* d (v^ dx 1 + v^ dx^ + dx^) 


= * 


3V 1 . 2 1 3v l 3 A . 1 3v 2 . 1 '2 

dx a dx + —— x- dx a dx + — =- dx a dx 

,3x .• dx J 3x 


Equation (A.12-3) continued 
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(A. 12-3) 

A. 13 Dot and Cross Product Identities 

Some of the identities involving dot and cross products will be 
listed here with proofs as needed. Since the dot product is defined 
by a symmetric bilinear form, 


A • B - B . A 


(A, 13-1) 


Similarly, since the cross product is defined by the exterior product, 
A*B=-Bxa. (A. 13-2) 

Using the tilde symbol introduced in (1,5—40), the cross product may 
also be written 

A x B = b Y e 

— — by a 


= [ e g“ e a® e <$] ‘[ bY e Y ] 


= A . B . 


(A. 13-3) 


Similarly, 


A*B = A.B = A*B = -B-A = -B.A=-Bxa . (A. 13-4) 

The vector triple product may be rewritten by applying Equation 
(A. 13-3) as 
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(A. 13-5) 


A • (B x £ ) = A • l' . ( A x b ) . C . 

The twelve permutations of (A, 15 , £) may be obtained from the symmetry 
properties, Equations (A. 13-1 and -2) and from (A. 13-5). All even 
permutations are equal, and all odd permutations are equal in magnitude 
but of opposite sign. In terms of the tilde notation, 

A.B X C=A*B.C=A.B-C=A*B- g , (A. 13-6) 

where the parentheses have been dropped since there is only one 
interpretation possible. 

The parentheses must be retained in the iterated vector product 
since there are two different forms. The iterated vector product may 
be written 

w /t, .. .a S . v .A r. 

(A. 13-7) 


A x (B x c) - e gy A-L . 


Noting the expression for the ‘Contraction of Levi-Civita epsilons 
(A. 9-14), this may be written as 

A MB x C) = (6“ 6 Se - 6“ ) a e b « ^ 


= (a " 6 3, cE)b V ^ « 6 6 » 6 )c a e a 

* (A . C)B - (A * B)C . 

From the anti-symmetry of the cross product 

A x (B x c) = (C x B) x A . 

Finally it may be noted that 


(A. 13-8) 


(A. 13-9) 


A x (B x c) + B x (C x A) + C x (A x B) 

= (A * C)B - (A.B)C + (B.A)C - (B*C)A + (C*B) A - (C-A)B 


= 0 


(A. 13-10) 
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In terms of the tilde notation, this last result may be written 


(g • B)~ • A = - B • (g • A) + Q - (B * A) , (A. 13-11) 

where the tilde outside the parentheses implies the quantity inside is 
to be operated on. This may be written in an operational form as 

/v * s “ l 

(C B) = C B - B C . (A. 13-12) 

The four fold iterated vector product may be expressed in either 
of two forms, 

A x [B x (C x D ) ] = A x [(B ■ D)C - (B • C)D] 

= ' D)A x £ - (B * C)A x D (A. 13-13) 

or as 

Ax (B x [C x £]) = (A • C x D)B - (A • B) C x d . (A. 13-14) 

A. 14 The Variation of a Vector 

The use of the variational operator, S, in conjunction with vector 

notation, in an inertial reference frame, has the interpretation 

6A = (6a“) - — . (A. 14-1) 

3x a 

This may be expressed in terms of the A ^ related basis, (1.5-35), as 


6A = (6a 01 ) c g - — ^ 
“ 3x' e 


■ *“ TP 


, a* 3> 3 
(a Sc ) s- 

a 3x ,p 


(A. 14-2) 


Ct 6 6 v "Y ft 

If the combination, a c , is denoted 1 a p (noting that a' = c* 1 f a p = 

a p 

3 3 

c^ ’a ), then Equation (A. 14-2) becomes 
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(A. 14-3) 


<5A 


(6 'a 3 ) ^ 


Sx 


-s- - ('a Y c Y 6 c e ) - 

,3 a a 


3x' 


» 6' A - ('a Y c Y 5 c 3 ) - — -r , 
“ “ 3x’ 8 


where the symbol 6 ' Is defined by 
<5 ' B = ((Sb a ) . 


(A. 14-4) 


dx' 


Y 8 

It may be noted that the array, ^ 6 , is skew symmetric since 


0 = 6(6 y3 ) = <S(c Y c*) = 5 c Y c 3 + c Y 6c^ 

a a a a a a 


0 


6c Y - - c Y 6c^ 


T a a a a 

Thus in a dot product, 


(A. 14-5) 


6 (A * B) = 6a a b a + a a 5b a 


= (6‘a^r- 1 a Y c Y <$c^) *b^+ 'a^d'b^- *b Y c Y 6 C e ) 


a a 


a a 


= 6 'a^'b^ 'a^S'b 6 - ' a Y b f ^c Y 6c^- *a P, b Y c Y 6c e 

a a a a 

= 6 'a'^’b 6 + f a 3 6’b e 
= 6 f (A * B) . 


(A. 14-6) 


In the case of mixed bases, 


‘a 9 tJ-3 1 «( a *3 , 8. 

r 17 j^J “■ '« “ > 


_ - a B Tl 3 , ou 0,,B , a 8 CI ,8 
- oa c T b + a 6 c *b + a c 6 *b 
ot a a 

= 6A • B -r a a, b p, 6c^ H A • 6 f B . (A- 14-7) 

— a — — ■ 
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A. 15 Modes 


An ordinary, second order constant coefficient system of 
differential equations in the form 

Mq + + G c[ + (1 + oti)Kq = 0 (A. 15-1) 

must be converted to first order before modes can be calculated. The 
matrices M, C, K, are real positive semi-definite, G is real skew- 

symmetric, and i = \/-l. First the q must be transformed into those 

1 2 

spanning M, q, arid those spanning the null space of M, q. Pre and 

post multiplying by the transform gives 


'll • 


1.. 


li , 

12 " 


1. 


ii„ 

12 ~ 


1. 


11, 

12 1 


1 

MO 


_i 


c 

C 


q 


G 

G 


q 


K 

K 





2.. 

+ 

21 0 , 

22 0 


2. 

■+ 

21„ 

2 


2. 

+ (1+ai) 

21, 

22 


2 

o 1 
o 

— 1 


q 


C ‘ 

C 


qj 


G 

G 


q 


K 

K 


q) 


(A. 15-2) 


The equation may be written in first order form as 


* li ' 1 

0 , - M » 0 


ri*o 

q 


13 1{' 0 . 0 


fl.l 

_q 


■< 

0 

11 M i 11„ 11 *12 12„ 

M G+ C i G+ C 

H 

l 

q 

►+ 

* 11 ’ 12 
0 (1+ai) K^l+ai) K 


l 

q 

K = « 

0 - 

o' 21 » 21 c' 22 » 22 c 


2 

q 


1 21 i 22 

0 ,(l+ai) K (1+ai) K 


2 

q 


0 

i 






- * 




(A. 15-3) 


If C is zero, then the first matrix is skew-symmetric. If G is zero 
(and C is not), then reversing the sign on the first equation would 
result in both coefficients being symmetric. 

A first order differential equation of the form 

A q + B q = 0 (A. 15-4) 
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may be subjected to further reduction if A is of less than full rank. 


Transforming q into those spanning A, q, and those spanning A's null 
2 

space, q, and then pre- and post multiplying by the transform gives 


n A' 0 


0 ' 0 



11 

1 12 " 


r > 

1 

r 

B 

B 

i 


q 

< 

21 b 

» ~ 

22 

i 


2 

q 

( 


(A. 15-5) 


The last of these equations is algebraic and can be solved for the 

21 22 22 

space spanning the partition B B. If B is non-singular this is 


simply 


-1 


N - - N" N H • 


which when substituted in give 


(A. 15-6) 


[ U a] { 1 q} + [ 1X B - 12 B ( 22 B) 1 21 b] j^-qj = 0 . 


(A. 15-7) 


The singular case yields similar results. In either case, we are 
returned to an equation of the form 

D4 + Eq = 0 (A. 15-8) 

where D has full rank. 

The procedure for reducing the dimension of the equations 
generalizes a method applied in structural analysis known as static 
condensation or Guyan reduction. Those methods were developed for the 
special case where moment of inertia terms are considered negligible 
and dropped from the mass matrix, reducing its rank to about half. 

In those applications the mass matrix is frequently diagonal, there 
is usually no damping or coriolis terms, and the resulting K is 
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non-singular , allowing immediate identification of the null space and 
solution in the form of Equation (A. 15-6). In the more general pro- 
cedure, finding a null space is not trivial. Eigen routines may be 
used (if they can treat multiple zero roots), some Choleski decompo- 
sition routines, based on maximal pivots, provide the spaces, and 
there are a class of elimination methods sometimes known as "structure 


cutters" which do the job. 

Returning to the Equation (A. 15-8), if E has less than full rank. 


transforming gives 


11 ' 12_ 1 11 * 1 
D , D q E 0 q 

1 i — r + — ' - v - 

21 n 22 2 ' n 2 
D i D q 0,0 q 


(A. 15-9) 


and solving as in Equation (A. 15-6) provides the rigid body modes and 


leaves 


Fq + Hq = 0 . 


(A. 15-10) 


The matrix H is symmetric and may be factored into lower, L, and 
T 

upper, L , triangular matrices (which is a by-product of- using a 

-1 -IT 

Choleski routine on E) . Pre- and post multiplying by L and L 


provides 


Jq + Iq = 0 , 


(A. 15-11) 


—l —It 

where I is an identity, and J = L F L . Assuming an exponential 
solution q = q e provides the standard eigenvalue problem on the 
dynamic matrix J 


+ A[Ij 


l] (q> = 


{ 0 } . 


(A. 15-12) 


* 
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The frequencies of the system are 1/A so the usual procedure of 
calculating highest eigenvalues first provides the desired lowest 
frequencies first (or, with a shift, in some low range). 

It may be noted that in the absence of damping (C = 0, a = 0), 
i[J] is hermitian, so the eigenvalues, iA, are real and the resulting 
motion is sinusoidal. In the general case, the adjoint problem must 
also be solved 

[[J] T +X [I]J { a 5> = Co} . ■ (A. 15-13) 

The eigenvalues are the same, and occur in conjugate pairs or are 
real, but the eigenvectors are different. After transforming the 
vectors back to the original q’s, premultiplying by the transpose 
of the adjoint vectors and substituting the vectors times modal 
amplitudes, [<J>] (nl, into the first order equations (e.g., A. 15-3) 
diagonalizes the system of equations. 
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APPENDIX B 


STRESS 

The presentation of the relationship between Newton’s Laws 
and stress has been removed to this appendix to preserve the 
continuity of the text. 
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B.I Stress Transformations 


Let m be an arbitrary point in a material continuum, N an 
arbitrary spatially contravariant vector of unit length. Let 


(x , x , x ) be Cartesian coordinates with origin at m, and such that 
the point defined by eN lies in the first quadrant. The plane normal 
to N through eN, and the coordinate planes define a tetrahedron with 



The surface force densities on the planes normal to N and the coordi- 

ct N ex 

nate axes x will be denoted jS and - _S respectively. The volume 
of the tetrahedron is 



(B.l-l) 

(B.l-2) 
(B. 1-3) 
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From Section 1.6, Newton T s laws imply 


= /• *(- ' i?) *L'k ' b) ■ 

Since the terms P, A , X remain constant as £ varies ^0(1)^ and the 
integral varies as the cube ^C>(£ 3 )) , and noting 


S_(m f ) = _S(m) + 0(e) , 
the equation becomes 


-“si-^gL + N S 1 e 


“2 12 3 

n n n 


-2 12 3 

n n n 


- 0(e . 


(B.l-5) 


(B.l-6) 


Thus 


^(m*) = n a °jS(m) + 0(e) y 

In the limit as e 0 this gives 

N c , v a a_, . 

_S(m; ■ n S(m) . 

B.2 Generalized Stress Transformations 


(B. 1-7) 


(B. 1-8) 


As in Appendix B.l, let m be an arbitrary point in a material 

continuum, N an arbitrary spatially contravariant vector with unit 
12 3 

length. Let (x ,x ,x ) be an arbitrary coordinate system such that N 
lies in the first quadrant. The plane normal to N and at distance £ 
from m, and the planes defined by the tangents to the coordinate curves 
at m describe a tetrahedron. 
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a/ax 1 


The vectors from m to the apexes may be denoted 

i a 


u = ir 


3x' 


(No Sum) . 


(B. 2-1) 


The product of the area of the triangle and a unit outward normal 

vector is given by the cross products 

1. 1 3 T 2 

A = ~2 U x £ 

2 2 T 1 1 3„ 

A N = y, U x U 


3. 3 1 2 1 

A N = Y U x 1 . 


For the face normal to N ( 


-A I = i( 2 U - \) * ( 3 u - A) 

= f( 2 £ x 3 £ + 3 £ x \ + X u x 2 u) . 

As in Appendix B.l, in the limit as e-»* 0 

N A N S e + *A 1 S^ = 0 . 

To demonstrate that this holds for 

N s = H ■ I 

Equations (B.2-2 and —3) may be substituted yielding 


(B. 2-2) 


(B. 2-3) 


(B. 2-4) 


(B.2-5) 
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ro|h-* 


(B.2-6) 


( 2 u x \ \ X ^ + ^ X 2 u) * £ + y| 3 U X 2 u 

"U x \ • £ j = 0 , 


. z 


, 1 1T 3 tt r , 2. 
+ U x U • jZ. ■+■ 


which is satisfied identically. 
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APPENDIX C 


BEAMS 

The calculation of numerical results for a beam finite 
element are presented in this appendix. This material is placed 
in this appendix to preserve continuity in the text. 


preceding 




C.l Beam Finite Element 


The material continuum is a beam of length depth h, and width 



The position vector will be approximated by 

P = [X + u(x, y)] ±-+ [y + v(x>] + [7] — , 

where the displacements are 


u 


V 


1 -2 


u + 10 U x + ni u y + T1 U x y + „ u - * y 


01 


11 


21 2 


0 v + lV X + 2 v i X 2 + 3 v i X 3 


(c.l-l) 


Surface forces will be prescribed as zero on the lateral surfaces of 
the beam. Body forces will also be prescribed as zero. The prescribed 
positions will be expressed in terms of nodes located at the center of 
each end. They will be assumed to have the form 


+ h ■ 


F -h=0 = (l q - 5 q y ) h + (y + 3 q ) !y 

P ilx=Jl = [ l + ( 2 q - 6 q *)] h + (y + 4 q )l7 + (7) h • (ca - 2) 

The element will be restricted to infinitessimal (or at least 
su ^ iciently small) strains so that the strain energy density can be 
approximated by a quadratic in the strains 

A - A + A E ag +- A E . . CC.1-3) 
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The material is required to be homogeneous, so that the coefficients 
are independent of the coordinates. And, the material is required to 
be isotropic so that the coefficients are determined by Young r s modulus, 
E, and Poisson's ratio, V 


Aim 

A 


1 + V 


. , AIIJJ , 

+ A , A = A , 


AIJJI AIJIJ , v /r i /\ 

A = A = G (No Sums) , (C.l-4) 

where I ^ J, all other coefficients are zero, and X and G are Lame's 

constant and the shear modulus defined by 

^ VE „ E 

A “ (1 + v)(l - 2v) * and G 2(1 + v) 

Taking the variation indicated in Equation (1.8-26) gives 


*h' 

-/<[('■ - ““] ■ 


The coefficient of 5e 0 

ap 


must be zero. 


so 


£ _ ^a£y6 £ 

y6 


(C.l-5) 


(C.l-6) 


The beam will be taken to be in a state of uniaxial tension along the 
— Axx A 

x-axis in the reference state. Thus A = a, and all others are zero. 
The term "VP" may be calculated from Equation (1.8-19) as 
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(C.l-7) 


( VP ) = 1 “ + i f / 3“\ 2 + /^\ 2 + / 1«\ 2 1 

-- 11 3 x 2 [Ux) w) IdJ 

/nr»\ l(~8u 3v 3u Su 3v 9v 9w 9w| 

(IE) 9 “ 71“ + ™ + + + 

[9y 9x 9x 9y 9x 9y 9x 9yJ 

2 > 2"1 

\ + /3w\ 

\9y' J 


*•»-$** KS *(| 

Inserting Equation (C.l-1) yields 


(VP) 


u +ii u y + ti ux y + 


11 10 11 ' 21 


t' 


+ (,v + 2 vx + 


2 [( 1qu 

1 — 2\2*1 

V 2 X )J 


— \2 


+ ll Uy + 21 UXy 


) : 


(VP) =- 
- L 'l2 2 


(VP) 


22 


[oi U + ll ux + 21 u 1 x ) + (l v V x V I x2 ) 
+ (l 0 u + 11 U V + 21 u X y ) (oi U 

" 2[(oi u + ll u x + 21 U 2 x2 ) ] 


+ ll u x + 21 l 


a3 


The coefficient of 6E p i n Equation (C.l-5) must be zero, 
only linear terms, this gives 
e. 


11 = 10 u + ll u y + 21 u X y 


12 


■ i[oi u 


+ 11 U x + 21 U I x2 + l v + 2 V x + 3 

Inserting this result in Equation (C.l-6) yields 


v I x2 ] 


A 11 

E 


22 


£ 33 - A 


a 12 

E 


= S+ (ttV +x )(io u 

= X {l0 U + ll u y + 21 U x y ) 
(lO U + ll u y + 21 u x y ) 

G (oi u 


+ 11 U y + 21 


u x yj 


+ HU x + 21 u j x 2 + x v + 2 v x + 3 v j x 



(C. 1-8) 
Retaining 

(C. 1-9) 


! j (C.l-10) 
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The third integral in Equation (C.l-5) may now be written as 


fs |[ S + (l^ + A ) (io u + n uy + 2i u 7 y)] [ 6 io u + s n u y + 6 2i u 
+ (io u + u u y + 2i u * y) ( 6 io u + 6 u u y + 6 2i u * y ) 

+ (l v + 2 v 7 + 3 V 2 * ) 

2 ° 
( 6 i V + 6 '2 V ' : + 6 3 V 2 * )] + [ G (oi U + 11 U 

+ 2 vx + 3 V | x )] [& Q1 u + « U u x + 6 n u \ x + 6 v + 


x y 


X + 21 U I x + 1 V 


6 2 v x 


, , 1 - 
+ 6 3 V 2 X 


f 10 U + 6 U u y + 6 21 U X y ) ( 0 1 U 


+ U u x + 


+ 6 1;L u x + 6 21 


( 6 01 U 

(oi U + 11 U X + 21 U 1 * ) ( 6 01 U 


1 - 2 \ 

U 2 X I 

(lO U + ll u y + 21 U X y ) 

u i x )] + [ x (io u + n u y + 2i u * y)] 


+ 6 n ux + 6 


u t x j 


(C.l-11) 


Collecting coefficients of the independent variations (i.e., 6 u, 

^ Ol 11 ’ ^ ji u » ^ 21 u * ^ l v * ^ 2 V * ^ 3 V ^ » retaining only linear terms, 
and defining 

f&y&z. = k,fy dydz = 0, J y 2 d^dz" = I, dydz" = h 2 I 


fdx = ^./x dx = 4 >/ x 


*2 - £' 
dx = 


Equation (C.l-11) becomes 
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[ d A£ + (l^ + X ) A * 10 u + ° A * 10 u ] 6 10 u 

+ [g a(j, 01 u 

+ [(iw + x ) 1 (hi u + Y~ 2i u ) + ° 1 ( £ n u + T 21 U ) + G A (l- 01 u + T n u 

4 


+ — u u + r 2l u + V + 7" 2 V + 6 3 V ) 5 01 U 


n. 4 it 2 £ 3 l 

+ r 2 i u + r i v + r z v + r 3 V 


+ 


i— 

\l+v 


)] 6 u u 

) I (|- ll u + f 2 l U ) + ° 11 U + I" 21 U ) + G A ( 


3 21 

£ 5 

+ 8~ ll u + 20 21 U + t~ l v + IT 2 V + 20 3 


r_ 

6 01 l 


A ( £ 1 V + 2 2 V + 6 3 V ) + G A (^01 u 


6 2l“ 


9 } 9 ? 

+ rn u + r 2 i u + V 


[5 # 1 


2 3 

.xr . r 

+ 2 2 V + 6 3 V 


V 


S, 3 2 4 

V + 3" 2 V + 8" 3 


•) 

v)+ G 




u + 


u + 


2 01 3 11 8 21 


i 2 i 3 i> 4 

+ ~ i v + r 2 V + t 3 V 


1 

d A (f 1 V + f 2 V + H 3 v ) +GA ( 

+ * 3 + * 4 + 2 5 ) 

+ ~i v + r 2 v + ™ 37 


S 2 V 


* 3 u + + * 5 
7 — U — 11 4- ■ — ii 

6 01 8 11 20 21 


V 


(C . 1-12) 


The two remaining integrals in Equation (C.l-5) are over the 
deformed (physical space) ends of the beams. A rectangular area 
perpendicular to the x-axis, with sides dy, dz in material space, 
becomes, as a result of deformation, a parallelogram with sides 

(1 + E^dy ^2* ^ + E 3^ dz —3* where 1 + E 2 « J 1 + 2^22 > and 
is a unit vector tangent to the deformed x^axis. The ratio of 
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the area after deformation to the area before is (1+E 2 > (1+E ) sin 

*^2 *— 3 ^ * w ^ ere cos (i 2 > — 3 ^ ~ £ 23^ * The Jacobian of the 

transformation is equal to the ratio of the areas 


J = <1+E 2 )(1+E 3 ) / 1 - 


'23 


(1+E 2 ) 2 (1+E 3 ) 2 




^ 1 + 2 £ 22^ 1 + 2 e 33^ ~ £ 23 (C.l-13) 

Noting Equation (C.l-9), to linear order J — 1 . Transforming the 
fourth integral in Equation (C. 1-5) to material coordinates gives 
(| 


fzs * p(V + 6 io u x 


+ 6 Q1 uy + 6 ljU x y + 0 21 


1 -2 -\ 

u i x y ) 


+ S y |( 


6 n v + (S.v x + A v ^ x 2 + 6 


O’ ' “ 1 " A T “ 2 * 2 x + ° 3 V \ x3 )] 


[J] 


(C.l-14) 


Evaluating this integral at the two ends (i.e. x = 0 and x = £) and 
introducing the symbols 

/sl x . 0 *( sX t j n - x q 
is| x _ 0 * lsX y [Jl} = - 5 q 
M*. 0 * {S5,[J1} = 3 * 

^! x= / {sX[J]} ■ 2 ^ 

[j]} 


- 6^ 


/as, *{s y [j]} = q 

JdS x=£, 4 


(C.l-15) 


yields for Equation (C.l-14) 
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( 1 Q 5 0 u - 5 q ^ + 3 q V> + ( 2 Q V + 2 Q M 10 u - 6 Q 6 01 u 

2 2 
‘ 6 Q * 6 11 u - 6 Q T 6 21 u + 4 Q V + 4 Q W 1 v + 4 Q T V 


+ 4 Q t" 5 3 v) 


(C.l-16) 


Summing Equation (C.l-12) and the negative of (C.l-16), and setting 
the coefficients of the independent variations equal to zero provides 
the stationary conditions 


V : - ^ - 2 Q = 0 

5 1() u : a AZ + & AZ 


U + (ll? + X ) 


1 +v + V M 10 u - 2 Q 1 ■ 0 


+ T 3 V ) + s Q V " 0 


/ J X~ 9T * 9° \ 

6 01 u = GA \*-01 U + r"ll U + 6"21 U + V + 2 2 V + T 3 V ) 

6 n u : Mhi u + l"2i u ) + (i^r + x ) 4 i i u + f 2 i u ) 

+ AG (t-01 U + Fll U + l“21 u + f~l v + f~2 v + ^3 V ) + 6 Q 1=0 

(foi u + I" 11“ + fo 2 1 U + T~ l v + 1 “ 2 V + fo 3 v ) 


+ AG 


+ 6^ T m 0 


V 


6 1 V 


: - 3 Q - 4 q “ 0 


■ Mv -f i- 3 -) 

» a (t 01 » ♦ £- U u 


3 2 

. il . 0 A A 

6 21 U + £ 1 V + 2 2 


2 3 \ 

v + r 3 V /- 4^=° 


Equation (C.l-17) continued 
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V : ° A (f l V + f 2 V + F3 v ) 

+ AG (r oi u + T 11 U + T 21 U + T i v + T 2 v+ r 3 v )-4 q I - = 0 

, A / 2 3 a 4 a 5 \ 

V = ° A \^ 1 V + S' 2 V + 20 37 

f 


, .Ji? . i A 

+ AG \6 01 u + 8 ll u 


5 3 4 5 

£ , £ J ^ £ , a . _ n 

+ 20 21 U + 6 1 V + 8 2 V + 20 3 V i“4 Q 6 0 


(C . 1-17) 

These nine relations may be solved for the six Q’s and three constraints. 

iQ = 


- P + aS 10 U + A fe + x ) io u ] 

,Q = A&+A$ 1Q u + A (-^+ x)i 0 u 

[a a L £ 2 ^ l 3 1 . AG Tn 

^ = - (l a [5l V + Tl V + r 3 V J + ~ L 01 U 


+ 2~ 11 U + IT 21 U 


0^ ■ n 3 3 

+u i v + r 2 V + r 3 V 


o = A 
4' £ 


% 


V + 2 2 V + 


2 3 

. I . A 

+ “ n u + r 21 U 


1“ 3 V ) + F|v 

+ (h V+ |- 2 v + T 3 V )] 

Q = f t ll u + f 21 U ) + I (lf^ + X ) ( 1 11 U + T 21 u ) 

a.g(i z , sf^_ . i?_ . ii „ . j£ \ 

' i \2 01 + 6 ll u + 24 21 + 2 1 + 6 2 + 24 3 / 

la/ ) I / - X 7 " 2 

| Q - JL ni U 2 2l“/ 2 i 


l+v 


AG 

£ 


(s? 

\2 01 U 


2 3 4 2 

* , r ^ £ 

+ ■■■■■ u 4* — 

3 li 8 21 


) ( X 11 U + T 2l u ) 

u + T 1 V + 3“ 2 v + r 3 v ) 


Equation (C.l-18) continued 


with constraints 
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^ 3£ 

2 11 20 21 U 


5 £ 21 U + (lk + A ) £ 21 U = AG (oi U 

i 3£ 2 ) 

+ 1 V + 2 2 V + 20" 37 

§ a( 2 v + | 3 v)+ A G ( n u + f 21 u + 2 v + | 3 v)= 0 
0 A (2 V + 15 £ 3 V ) + AG (ll U + 15 ^21 U + 2 V + 15 £ 3 V ) = 


the last two of these imply 


§ 2 V + g( u u + 2 v)= 0. or 2 v = - U u 
$ 3 v + G ( 21 u + 3 v)= 0 


a+G 

G 


or 3 V = “ A 

a+G 21 

and inserting these in the first of Equations (C.l-18) gives 

aSg £ 


v i = - oi u + b - 


AOG 3£ 2 , a , 

2 11 U A, „ 20 21 U + 9 X 21 U 


a+G 


(l+V + X ) X 2l u 


+ \i+7 + x 

a £ a' 3£“ 

u - — - ^ u - 


a I 

u + — — - „ , u 


01 S+G 2 11 S+G 20 21 G A 21 


+ /2U-v]\ I 

+ \[l-2v]J A 21 U 

Substituting back in Equation (C.l-18) yields 

! Q = - A ° - A ^10 U - A (l^ + A )l0 u 

A$ 1Q u + a(— + xj 1Q u 


,Q = Aa + 


Equation (C.l-21) 


(C.l-18) 


(C. 1-19) 


(C. 1-20) 


continued 
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A = 


& 2 


o $ 2 to 2 

Ao' u + A — _ u + A - — -j—- ^ 

01 §+G 2 11 S+G 20 21 


A2 t 
. a I 

U - A G- I 21 U 


A a 


2 (1-v) 


(l-2v) 


I 

A 21 U 


, A CG_ i 

a+G 211 $ + G 20 21 


AaG 3£ 2 a t / E , 

U - &I 21 U - + X 


) l 21 U 


.Q = 


01 


&+G 2 ir 


- AoL. u - A - — ~ ,,u - A 


A2 oo 2 A2 t 

a 3& . . a I 

a +g 20 21 u A G A 21 u 


+ 


. A ^a-v)] i 
Aa Ll-2v J A 2l u 


A<JG l 

A&G 

9 

3 Jr 

A 2 

11 U A 

20 

O+G 

a+G 



„ 2 \ 


I<* la 

t £ \ 

I 

5 Q z \ h u 

+ t 2i u ; 

+ i 

, A „ 2 

/\ ■ 

3 

AaG it 

, AaG 

3-i 

\ ^ 

11 U + A 

40 

a+G 

O+G 


. A ,2 

A 

3 

AaG i 

AaG 

11 


h + A )( a U u + f- 21 u ) 

U - Sl f 21 u - 21 u 


•Q = 


$ + G 6 11 &4G 24 21 


" T { l U U + 2 21 u ) ~ l {lh + A ) (*• ll u + f - 21 u ) 

+ fS T n u + rr Ir 2i u - Sl f 2i u - fe + A ) T f 21 “ 


a+G 

AaG J2 2 . 
§+G 3 X1 ‘ 


a+G 

aSg 

$+G 8 21 


(C. 1-21) 


It will be assumed that the preload is very small compared to 
the shear or elastic moduli (e.g., for steel G = 11000 ksi, E = 29000 
ksi, and the yield stress is 36 ksi; so a < .0033G to be in the elastic 
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range). It will also be assumed that the beam is slender (e.g. h << Z) 

2 


121 3 

so that • — *r « 1 (recall I = wh / 12, A = wh, 

Al 

(C.l-21) can then be written 


so 


121 
A Z 2 


-«)>• 


Equation 


= - Aa - Aa^u - AE 


? Q = Aa + Aa 1Q u + AE 


[(1+v) (1-2 v) ] 10 u 

L-v 


[(1+v) (l-2v) ]l0 u 

,2 


> Q= A %1 U + I 11 U + §~ 21° - EI [ (1 +v) (l-2v)] 21 U 


4 Q « - aS q 1 u* - Act f a - a 5 3 * 


u + EI 


I— l=v 

L(i+v)(j 


20 21“ ’ [(l+v)(l-2v) 

-v 



5 Q = + Ao — U u + Aa ^ n u + ei [(J^' ) ( 1 _ 2v) ] 1L u 

6 Q “ ‘ A0 12 ll u “ Aa 20 21 u " EI [u+v) V (l-2v)j ll u 

~ E1 [(1+V) (l-2v)J 4 21 U (C.l-22) 

The above noted approximations along with the constraints (C.l-19 
and -20) may be used to re-express the second of Equations (C.l-1) 

v = o v ~ oi u * - n u f f * - 21 u f ir * ■ n u t * 2 - 2i u £ * 3 

(C . 1-23) 

The final integral In Equation (C.l-5) when converted to material 
event space assumes the form 
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i 


3S| X=0 * (|(o U + 01 U y - ! q + y 5 q ) 6S + (o V - 3 q ) 

+ + ° U + lt)Ue + <>l Uy + ^ + 2lU ^ 

‘) 


«S y 


* 2 - 


- & - 2 q + y ^)& s ' 


+ lo V " 01 u£ “ 11 U G 2 21 u G 20 “ ll u 2 


5 3* 3 




u *£“ “ A q)6S y |j| 


21 6 4 

The stationary conditions are 


(C.l-24) 


SS - n : 
1 x=0 

0 U “ l q = 0 


x_ 

Ss y lx-0 : 

01 u + 5 q = 0 


6S y |“ . : 
1 x=0 

0 V - 3 q " 0 


6S X |- . : 

1 X-Jt 

0 U + 10 Uil - 2 q = 0 


6sXy lx=Jl : 

0l“ + ll 11 *' + 21 U 2 + 

6 q = 0 

* 

j, 2 

0 v " 01 uA 11 U 2 

i 3 

21 U 6 


q = 0 


(C.l-25) 


Solving for the displacement coefficients 


o u = i q 


10 


u - - I l q + I 2 q 


01 U = - 5 q 

6 6 ,4 2 

H U = ^2 3 q ~ ~2 4 q + ~ 5 q + I 6 q 

, _ 12 12 6 6 

21 ' ' * 3 3<J . * 3 4<1 ' t 2 5<1 ' ? 6 q 


0 V = 3 q 


(C.l-26) 
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Substituting this into Equation (C.l-22) the result may be expressed 


in matrix form as 


where 


K g = AO 


Q = Q 0 + K G q + Kq 


Q - [jQ, 2 q, 3 q, 4 q. 5 q, 6 q ] t 
Q - [ - A a, aS, 0, 0, 0, 0 ] T 

T 

9 “ t » 2 q * 3 q> 4 q > 5^* ^9 1 
0 ! 0 


(C. 1-27) 


K 


El 


^ 1 

f 

~ T ~ ” 


3/A - 9/51 l ~ 3/A + 9/5A 

I 

= 6/5A [ = - 6/5A 


- 3/A + 9/5A , 3/A - 9/5A 


= - 6/5A , 


= 6/5A 


1/2 - 2/5 1/2 + 2/5 

i 

- 1/10 1 - - 1/10 
- - - - - + - - - - - - 


- 1/2 + 3/5 
= 1/10 


1/2 - 3/5 
= - 1/10 


- 1+2 - .9 1 1 - .9 


= 1/10 = 1/10 


1-2 + .9 f -l + .9 
= - 1/10 \ = - 1/10 


A/3 - A/5 1 A /6 - A/5 

\ 

= 2A/15 f = - A/30 

- - _ _ _ _r _ 


-A/3 + 3A/10| -A /6 + 3A/10 

= - A/30 1 = 2A/15 

f 


AA , AA , 

i r i i 

0 | 

0 1 
i 

0 

1 0 
i 

AA 1 AA 1 

"T 1 l“' 

o ! 

\ 

0 , 

0 

i 

1 0 

0 i 0 i 

12/A i 

- 12 / l 1 

6 

i 6 

~o~ 0 2 ~ 

" - 12?A“ T ~ ~ 

~ ^?/r ' 

~6 

- -_- 6 

0 i 0 i 

6 "I 

- 6 i 

4A 

i 2A 

0 1 o 1 

- 6 + 12 1 

i 

6 - 12 1 

i 

- 4A + 

6A 1 - 2A + 6A 

1 i 

1 1 

- 6 1 

— I 

- - 6 1 

= 21 

; - hi 


where E = E 


( (1 + 


1 - V 
V)(l - 2v) 


)• 
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This result is identical to the geometrically nonlinear beam 

finite element usually obtained by assuming a stress distribution, 

except for E replacing E. This is due to the fact that the simplified 

displacement field used in this example allows no extension in the 

lateral direction, and the associated stresses fail to vanish on the 

lateral surfaces. This deficiency may be corrected by selecting a 

more complete displacement field, or merely by replacing E by E, 

noting that E is the axial modulus for laterally restrained material, 

and E for laterally free material. That is, E is the elastic coefficient 
XIII 

A (C.l-4) giving the stress for unit strain with the other normal 

Jill 1 


strains zero. The inverse set of elastic constants B' 


E ’ 


= - — provide the strains given the stresses. Thus for a unit 
h 

stress with the other normal stresses zero, the strain is 77 , 

E 

C.2 Beam Dynamic Parameters 

The purpose of this section is to develope the dynamic parameters 
for the finite element described in the previous appendix (C.l). The 
assumed position field may be obtained by substituting the results of 
Equations (C. 1-19 ,-20,-26) in Equation (C.l-1) 

P = 


X + l q + X H l’ + I 2 q ) + y f 5 q ) 
+ * y /• 


6 

4 , 2 \ 

3 q - ? 4 q + 

1 5 q + 1 6 q ; 

1 12 12 

6 6 \ 

* 3 3<5 £ 3 

4 q -^ 5 q --> 6 q ]j 


3_ 

3x 


— , , — / \ .1 — 2 / 6 6 4 2 \ 

y + 3 q + X ( 5 q) + j x — 3 q + _ 4 q - ¥ 5 q _ x 6 qj 


+ 6 X 


1-3 / 12 


(■ 


12 


33^ q ^ o n aQ 


3 4 


2 5 


2 6 


■) 


JL 

3y 


(C. 2-1) 
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From this, the functions £ qP ' , ^P* may be determined 



Assuming a uniform density p the dynamic parameters may be 
computed as 

M = * {p} (C. 2-3) 

= p w h £ 

The center of mass is 



The shift in mass center due to generalized displacements are 
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(C.2-4) 
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The moment of inertia is 




The change in moment of inertia due to generalized displacement is 




fs * | P [(eO^* ’ 2- ? ) H “ EO ^- 1 | 


Equation (C.2-6) continued 
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3 i =y^ *jp [( eo p’- 3 £')m - 

- aI u /ii i ^ ® (id -h 3 )/6 v 

- P [- Mh U - kl 2 5l 3j 3x ® 3y - W |24 ' 24 [(,2 2 

6 Z 3 \ 3 _ 8 1 

^3 3 / 3y 3xJ 


- m ! ii iL ^ 9 h2 8 ^ 3 < 

“ M } “ 20 3x ® 3y “ 121 3y ® 3x ' 


p* p* 
EO— 4— 


ii = /?*j[(Eo£’- 4 £')H- 

a[" X L ** 0 £ 5 \ 3 ~ 3 J h3 -h 3 !/ 6 ^ . 6 £ 3 \s ~3 | 

T ■" i 3 ^ - 2 - a-|(r rj» ®d 


5V 

J 7i3 . 3 x h 2 3 ~ 3 

' j 20 3x ° 3y + 12 i 3y ® 3x 


P ’ P ' 
E0“ 5— 


■]| 


->i = fs* {[(hoE’-^'Jh - 

2 £ 4 , 1 1 _ to 9 _ 

~ P l_" Wh \3 “ £ 4 + ^2 5 / 3x ° 3y 


__ M j £ 2 3 ~ 3 ) 
M | 30 3x®3yj 


fh 3 -h 3 
" W | 24 “ 24 



9 


3\ 

K— i 

% 

3_ 

e\ 

2 

l 2 

3 / 


3 ^ 3 

3 F® 37 


£ "fi * {[(ecT-^H - E W]} 

A . I n 4 ,l £ 5 \ 3 _ 3 

= Pl_ wh + 

Jhi -h 3 ) 1 2 l 2 3 £ 3 \ 3 „ 3 1 

W |24 ~W\ In 7 -) 37 ® 3 ^J 

= J + 

M j + 20 3x U 3yj 


(C . 2 - 6 ) 


The first order moment terms are 
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0 = /-*!p _P t x P»| 

1- Js r eo- i- j 
" 0 

2— ‘/s *|p EO— ' * 2— * J 
- 0 





It may be noted that these calculations are of the same form as some 

of the i I terms and need not be separately computed. The second order 
changes in moment of inertia are 
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Equation (C.2-8) continued 
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I =/s *je[( 2 r- 2 l , )“- 2 ^' 2 ^'] [ 
i =/s *{p[( 2 I’-3£’)m- 2 I’3i']} 


-- nt Uh®h 


24i =/?*{p[( 2 r*4l’)y- ^V]J 


M |lo| lx ®ly 


25i = /s «’]! 




33= 


S“/s *'p [( 2 £'-^ , ) 2 - £'£']J 
= + “| 20 1 
■/s*jp 


= Ml i.® ».+ ji-4 + ^! 

1 



Mil. <a JL £1 1 

35j8x ^ 9x | 10 z 2\ 9y ^ 9y T |l0 ^2 T 35j 9z ^ 97 


34= \A*jp ^i']) 


= M 


|2_| Il_ h^l 3_6a3_. ( i_ h* ^ 9 ) 3 „3 _ 

'70) 3x® 3x ‘ 1 10 £ 2| 3y ° 3y + 10 £ 2 + 70| 37 ® 37 


35= ■/? *{ P - a^'jJ 


= M£ 


11 ) 3_^3 , ) 1 h (3 ~ 3 
210) 9x U 9x + 120 0 2 9y® 9 y + )i20 


i_ h^_ u) 3_ n 1_] 

20 £ 2 210) 9z ^ 9z| 


Equation (C.2-8) continued 
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i=/ ? *j A p[( 3 £ , - 6 p , )y- 3^']) 

®| =+ !^4flr®lr + 


» Mil 


1 h 13 


!fe®sl 


'420) 3x 3x (120 l\ 3y U 3y (120 ^2 “ 420 

44 i =y^*{p[(^'-^’)2-^'^']j 

45±-y?*{ S {(4 , ‘5£ , )a-4£ , 5 £ , ]j 

- mo 1113 ) 3 ~ 3 ) 1 h 2 ) 3 ~ 3 . j 1 h 2 ,13 ) 3 -.3 

Mi |j420|3x ° 3x " |l20 ^2j 3y ° 3y + |~ 120 % 2 + 420J 3z®3z 

1 =7i- *j $ [(4-£'- & ^')= “ 


46 


- Mil 


3 , 

( 1 

h 2 , 

n ) 

3 ^3 

3y 

(120 

~ + 

l 

120) 

Cl ® a 

3z 3z 


53= -/? *js [( 5 f’* 5 £)a- ^’J 


(105) 3x ° 3x (90 £ 2j 3y U 3y j 90 ^2 + 105) 3z ° 3z 


i^*jp[( 5 P'- 6 P’)u- 5 P' 6 p| 

Jj 3_aL.l 

60 £ 2) 3y w 3y (360 ^2 l40)Bz U 3z 


-^ 2 |ilo!i®i + 360 


66 




p * P ’ 

6— (r- 


- MiT 


jJLj !_<*) 3 li-Jiil li_ h 2 , 1 h ~ 3 ] 

Il05) Bx ° 3x (90 ^2) By Q By + (90 ^2 + 105) 3z ® 3zj 


(C. 2-8) 

where the remaining terms may be determined by the symmetry relation: 


I = x a ^ i — <g) 1 — = 

ab= ab 3 x a 3 J ba 


X^ a 1 — 0 2 

. 3x a 3.x* 


(C. 2-9) 
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, - _ Mil 3^3 

(e.g., 35 J ' 53 J but 16 J 3Q 3x ® 3y and 61 J 


Mil, ® 1-), 

30 9y ^ 3x } ' 


Similarly, the second order moment terms are 


13— 


14- 


= M 


J A 

*jp 

P* 

1- 

|7J 

1_ 

20 ) 

3z 

* ( A 
* P 

P f 
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3 1 

3 

20 j 
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15— J S j P 1- 
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x p 1 

5— 


16° 


23 


- f— A j A P » X P * | 

J S | p 1- 6^- j 


“ - M y It 


o=y-*| P A 2 p< . 


-lloli 


24 


- 7 s * ! p -- x <■* 


« M 


j M 2- 
7 ) 0 




23 i. -/ s 


20j 9z 

* I n P f x p T 1 
( P 2- fr- j 


■ “ \M h 


26^ fs 


*jp 2 P- x gP-j 


= - m jloli 


(C, 


2 - 10 ) 
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The elements with interchanged subscripts are obtained from 

antisymmetry ( ,0 = - , 0). All other elements are zero. The mass 
ab — ba— 

matrix elements for a consistent formulation (i.e., assumed velocity 
field the same as the assumed displacement field) are 


, ,M = /— * J n ,P’ * P f [ 

= M |~| 

11 J S 1- 1- J 

(3 j 

12 M= /s 

= M |II 

13 M =/s *(P ll’-st’) 

= 0 

u M = fs *{" 1 £ ’ 1 

■ = 0 

15 M “/s *{P if'st’) 

= 0 

16 M = fs *| £ 

= 0 

22 M ” A *{ P 2—’ * 2—' ] 

• M &! 

23 M = fs * [ p 2-' ' 3—* ] 

f- 0 

24 M= /s *| P 2—' * 4—' 

= 0 

25 M =/s *| p 2£'-5t’) 

► = o 

26 J S [ P 2- 6r~ 

l-° 

33 M = /s *( P 

1_ M I 1 h2 

r |io 72 - 

34 M = /s *f £’*£’ 

}- M 

J | 1° £ 2 

M = /— * 0 P’ • P* 
35 y s ¥- 5- 

)■ » ji» 

36^ y S l p 3^ 6^ 



11 i 
210 ' 


420) 

Equation (C.2-11) continued 
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11 
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M 

55 

fi *C» 

p» 

5±_ 

5— * 

= m 2 

1 

9C 
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■fs‘ I s 

5^- ? 

•<£\ 

2 

= m 

-1 

360 

h 2 

£ 2 

1 

140 

66 M 

•fs'p 



2 

= MJT 

1 

90 

h 2 

£ 2 



where the remaining elements are determined by symmetry (i.e., 

“ j i^) ■ ^ may be noted that if the beam is slender (y « 1), then 

the terms involving (h/£) are negligible. The mass matrix is then 
identical to the consistent mass matrix for a slender beam as presented 
in the literature. 

C.3 Substructure Modes 

In this section modes for a very simple substructure are calculated 
using Hurty's method. The substructure selected for this example 
consists of a single beam finite element in a stationary reference frame 
(see Sections C.l and C.2 of this appendix). It will be assumed that 
the beam is attached to other substructures at ^q, 2 q, and 5 q. Out 

of plane deformations are not included and the preload is assumed to 
be zero. 

For this special case the equations of motion (2.3-5) become 

j* M $ + j£ K *9 ■ jQ (C. 3-1) 
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where M and K are obtained from Equations (C.l-17) and (C.2-11), 
respectively. Rearranging the rows and columns to 4,6,1, 2,3,5 to 
group internal (4,6) and boundary (1,2, 3,5) nodes, the equations may 
be partitioned 
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The modes of the internal degrees of freedom satisfy 

[ K II M II “ i* U ] {i^} = 0 (C. 3-5) 

2 

where =l/^cn , the inverse of the square of the frequency in radians 
per second, and ^<j) is the associated mode shape. This yields the 
characteristic equation 
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where A' = — r A. The roots are 
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with associated frequencies (radians/ second) 
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It may be noted the first frequency is quite accurate (the exact 
continuum mechanics solution is 3.52 
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Now assuming only the lower frequency is of interest, or merely 
noting that this single element model cannot accurately reflect 
behavior for higher frequencies, the higher mode is discarded. Thus 
the original coordinates are related to a new set of coordinates by 
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Transforming the mass matrix to these new coordinates gives 
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The reduced dimension mass and stiffness in an equation of the 
form (C.3-1) provide a reduced dimension set of differential equations 
in the modal amplitude, r\ s and the four boundary coordinates which 
describe the low frequency behavior. The reduction is not particularly 
impressive in this example since there are more boundary than internal 
nodes. However, in practice, the number of internal nodes is usually 
much greater than the number of boundary nodes. Discarding half of 
the modes in such cases results in substantial savings. 
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APPENDIX D 


QUADRATIC TERMS 

This appendix demonstrates the adequacy of a linear expansion. 
The unimportance of quadratic terms is demonstrated both analytically 
and by example. This material has been removed to this appendix in 
an effort to preserve the continuity of the text. 



D. 1 Quadratic Terms in Lagrangian 

The purpose of this section is to analytically demonstrate 
that in a Lagrangian formulation truncation of an assumed quadratic 
expansion of the displacement field to linear order results in 
equivalent (although not identical) linearized equations. Let the 
assumed displacement field be of the form 


P = + .P .q + . .P .q .q . (D.l-1) 

— 0- i— x ij— 3 » 

where the repeated Latin subscripts imply summation over all the 
generalized coordinates. Similarly, body and surface forces may be 
expressed in terms of the generalized coordinates as 

£ = qS + ± S ± q . (D.l-2) 


The generalized forces may be computed, 
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(D. 1-3) 


Selecting an inertial frame of reference (to avoid the detail of a 
quasi-coordinate formulation) the velocity and acceleration are 
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The kinetic energy may be calculated as 
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Lagrange’s equation (including all potential energy terms in the 
generalized force) , 
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Noting the cancellation of the third and tenth, fifth and eleventh, 
sixth and twelfth, and seventh and thirteenth terms, this becomes 
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Substituting the expansions in Equations (D.l - 1 and - 2) into 
Equations (D.l-3 and -8) and combining these two provides to linear 
order 
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Equation (D.l-10) continued 
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The zero order term is 
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Mow either by requiring the be a complete function space, noting 
the region is arbitrary, or merely recognizing the above as P times 
Newton’s law, it may be seen that 



must also be satisfied. Subtracting this from the linear term in 
Equation (D.l-10) gives 
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which is the equation obtained by retaining only linear terms. Thus 
the retention of quadratic terms adds the terms in Equation (D.l-12) 


which are redundant with the zero order terra (D.l-ll) and thus may be 
eliminated. 
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D.2 Quadratic Terms in Gimbal System 

The purpose of this section is to demonstrate, by example, 
that it is sufficient to retain only linear terms in an assumed 
quadratic expansion of the displacement field. To provide a reason- 
able test of the hypothesis, a system was selected for which rotations 
are the natural coordinates. The base of the illustrated gimbal 
system is experiencing a steady inertial spin, Q. The equilibrium 
state is 8^ « 0, in which steady state torques are carried through 


clocksprings. 




The linear equations governing small motions about the steady state 
equilibrium will.be obtained in three ways: first in the natural, 

Eulerian fashion, second retaining quadratic terms in a Lagrangian. 
formulation, and finally retaining only linear terms in a Lagrangian 
formulation. It will then be shown that the equations are all 
equivalent . 

This problem is not entirely academic. Structures are frequently 
represented by finite elements whose generalized displacements are 
three translations and three rotations at connecting nodes. Kinetic 
energy is sometimes associated with these models in the simplest 
possible way, by assuming all the mass is concentrated in inf initessi— 
mal rigid bodies at the nodes (which, however, also have rotary 
inertia) . The gimbal system may be thought of as a mechanization of 
such a nodal body. In this case the rotations are a 1-2-3 set about 
the body fixed axes. An alternative mechanization that has been 
suggested by other researchers is a 1-2-3 set about axes fixed in the 
substructure reference frame. It can be shown that such a set is 
equivalent to a 3-2-1 set about the body fixed axes, so the example 
covers both types of mechanizations. (Of course, the question of 
mechanization is academic, since it is about to be demonstrated that 
linear terms are sufficient, and since it is only the quadratic terms 
that are affected by mechanization.) 

Before proceeding with the development of equations, it will be 
convenient to express some of the relationships common to all of the 
derivations. It may be noted that the customary letter with a caret 
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(e.g., e) notation has been used for unit vectors in place of the more 
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(D. 2-1) 
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where the abbreviations c8 and s0 have been used for cos(0 ) and 

a a a \ 

sin (® a )» The angular velocity of the body is 


<0- 2 T n+-a, '6, + S* S, + b* 8, , 


(D. 2-2) 


-1 "1 ' ~ 2 v 2 ' ^ 3 w 3 
where the absence of subscripts in the first term implies the complete 
array treated as a matrix, and the juxtaposition of the two terms 
implies the matrix product. Making use of the direction cosine 
relationships (D.2-1) this may be written 
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In the body fixed basis this is 
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The angular acceleration is 

» At •• r* . »T — c— . 

OJ = V (0 + U 9) - a (0 + Q 6) . 
The components of the torque are 
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Expressing _T in the basis , 
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Multiplying the second equation by cO^ s0^ and adding c9 1 times the 
third equation gives 

ce 2 s6 l ( V K 2 6 2 ) + c9 1 (T 3" K 3 0 3 ) - ' c0 l se 2 (T r W = c6 2 T 3 * 

(D. 2-9) 

Linearizing this result provides 
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So the torque may be expressed as 



(D. 2-11) 

Before proceeding with the Eulerian derivation it will be 

demonstrated that a 1-2-3 rotation about axes fixed in the basis is 

the same as a 3-2-1 about body fixed axes. A set of three unit vectors, 

ju a , initially parallel to the a basis vectors, u - U a (where U is 
a a — — 

an identity matrix) , after a rotation 0^ about the axis have an 
expression in the a basis 

a .l = U 0(6^ a . (D. 2-12) 

A subsequent rotation, 0 2 , about the a_ 2 axis yields unit vectors, b ,fi, 
with an expression in the a basis 

fcf u =» U C(0 1 ) C(0 2 ) a . (D. 2-13) 

A final rotation, 0^, about the a_^ axis yields unit vectors, ^5, with 
an expression in the a basis 

Ji' u C(0 1 ) c(e 2 ) c(e 3 ) | . (D. 2 - 14 ) 

It may now be recognized that the product direction cosine matrix is 
the same as that for a 3-2-1 rotation about body fixed axes. 
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The equations of motion will now be developed following the 
usual Eulerian formulation. The rotational equation is 


T - *{p P x £} . 
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where I = I _b and I is the array of moments and products of 

A 

inertia. Writing this result in the b^ basis and noting Equations 
(D. 2-5,6, and -11) 
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Retaining only terms to linear order 
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These equations may also be written in the a basis. The inertia 
dyadic becomes 

X = a T (E +9) I (E - 0)a . 2 T (I +'01-10)2, (D.2-21) 

where the linear approximation to the direction cosine matrix (E - 0) 
has been used. The governing equations are 
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The steady state and deviational equations are 
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(where use has been made of the identity (A. 13-12)) and this is merely 
a restatement of the steady state torque relation. 

The equations will now be derived using Lagrange's equations 
retaining all quadratic terms. The notation is taken from Roberson 
and Likin's "The Quadratic Approximation in Rotational Dynamical 
Equations". With this notation the direction cosine matrix is 

A -J /W 1 TOO A 

b - (u - e + j e e + | e 6) a , (D. 2 - 26 ) 

where 
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The relative angular velocity, w, is expressed in the b basis as 
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The generalized forces may be computed as follows. Noting that 
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therefore 
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Thus the generalized forces are 
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and the potential energy is 
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which is just the sum of the strain energies in the three springs. 

The kinetic energy is 
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where it has been noted that x y = y x. So 
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Now noting x y = y x, the elements of Lagrange's equations may be 
obtained 
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|| - - Hi + y(u i23 ifi) A Je - Tim + Tim 

+ i S' (in>~ + | (in>~ n - - cTu 123 (iftr + 1 (inr u 123 n]e . 

(D.2-41) 

Lagrange’s equations may now be written 

t - kb « oio + iS +[1?? + Hi - (ift)~]e +[H i H 

- i {?T(ia)^ + + ~ {TTu 123 (iQ) " - (io) “ u 123 a}]0 . 

(D. 2-42) 

Thus the steady state and deviational equations are 

t = n i 

10 + [if? + Hi - (ift)~)e + fe 1 S + r - - {H(iq)~ + (i^nj 

+ i| 0 U 123 (If2)” - (IQ)“ U 123 H|Je * 0 . (D.2-43) 

Some of terms in the stiffness matrix can be identified by noting 

- ~ {H<ifo~ + (inrSi = - \ {2H(ior- mm)~- (mfT]} 

= - H(io)~ + ~(q m~ 

= - ff(Ifl)~ + ~ T , (D. 2-44) 

Subtracting the a basis Eulerian equations (D.2-23) yields 


222 



to J t— * 


0 - 

T 3 



0 

0 

■ 

0 

T 

3 

0 

~ T i 


- T 

3 

0 

0 

- T 

2 

T 

1 

0 


T 

2 

~ T 1 

0 


— 123 

[m J (iu) 


(IQ)' U 123 Q] 


0-0 


(D. 2-45) 


or 


0 

- T, 


^ (Ji 


0 - T, 


- 1 r )'X 101 

= (IQ) u ^ u - QU ^(IQ) 


T„ - T. 


(D. 2-46) 


Thus if this equation is satisfied Euler's equations and Lagrange's 
equations retaining quadratic terns are equivalent. That this equation 
is true, may be verified by direct expansion. From the steady state 
equation 


T = e ■ Q 0 I x Q. , 
a cxBy 3 y<5 6 

* T. « Q 0 I_ Q - Q I Q , 
1 2 3a a 3 2a a 


T 0 = Q 0 I, Q - Q. I_ Q , 
2 3 la a 1 3a a 


T a - Q_ I 0 Q - Q_ I- Q . 
3 1 2a a 2 la, a 


So 


(D. 2-47) 


o -t 3 

T 2 


0 

1 VlaV 

Q-I« Q 
1 2a a 

IQ-I- Q - Q-I_ Q 
( 3 la a 1 3a a 

-t 3 0 

- T 1 

“ 

Q_I, Q - Q.,I„ Q 
2 la a 1 2a a 

1 


1 Q I Q - Q.I. Q 
' 3 2a a 2 3a a 

1 

T 2 - T 1 

0 


Q - Q,I 0 Q 
3 la a 1 3a a 

, n o L Q - 
1 3 2a a 

Q 0 I~ Q 
2 3a a 

1 

i 0 

- 


(D.2-48) 
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The elements of the other terms are 
_ 123— 123 

{(xi2) ir = i aY 

^3 ^ no sum and 

{« U 123 (I«-} ae = tl a v£ 3 I YB n B (no sum a, 8). 




Thus they produce matrices 
0 


(IQ)"U 123 Q = 


a u 123 <ift) 


u 01 ^I 0 ft 
21 1 2a a 

u „ 1 ft.. I ~ ft 
31 1 3a a 


0 


ft 

21 2 la a 

u 01 ft.-I ft 
31 3 la a 


I u ft I ft 

12 2 la a 13 3 la a 

0 | u ftl 0 ft 

1 23 3 2a a 

0 


ft 

32 2 3a a 


^1 oft -1 1 a ft 
12 1 2a a 


u. 0 ft„I 0 ft 
32 3 2a a 


and 


U- ^ft 1 1 ft 
13 1 3a a 

U AAft 0 I„ 

23 2 3a a 


(D. 2-50) 

The difference between the two arrays, noting that U 123 - - U 323 eives 

ap 3a 


"" In ft ) 

12 2 la a 1 2a a 

Ui 9 (ft X ft — ft I ft ) 

13 3 la a 1 3a a 


U 12 ( Vlcftf ^lhaV ' “l^VlctV 


I u (ft I ft - ft I ft ) 
23 3 2a a 2 3a a 


23 3 2a a 2 3o a I 


(D.2-51) 

Recalling that u 12 = = u^ * 1, it can be seen that the arrays in 

Equations (D.2-48 and -51) are identical, so Equation (D.2-46) is 
satisfied and the Eulerian and quadratic Lagrangian equations are 
equivalent. It may be noted, however, that the stiffness matrix in 
the Lagrangian formulation is symmetric and its elements are order of 
rotation dependent. 
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The equations of motion will now be obtained using Lagranges 
equations, but retaining only linear terms in the displacement field. 
So in this case the position is 
P = C T P 

c£ - c o p 


AT «"•«' 

= a (E + 0) O P 


(D.2-52) 


The generalized forces are computed by noting that 

aT *** 

6 qP = a S6 Q P 


(D.2-53) 


SO 


Q 60 = f 

a a J a 


* 16 p - q) 
9M(t) 10 (£ 


69 { o p 


(D. 2-54) 


The velocity is 


(£-■£*■(£ + R x <<£ + I x <£) 


(D. 2-55) 


- (ft + 6) x p + x (0 x p) 

— — o— — — o— 

so the kinetic energy is 

T * *{pC (£ + 1) x qP + ft x (ex qP) ] • t (ft + 0) x Q? + £ x (£ XqP) ] } 

= *{p[(H + 9) x oP - (ft + i)x qP + 2(f2 + 0 ) x ^ . £ X (e_ Xq P) 

+ g X ( 0 X qP ) * £ X (Q X £ ) ] } 

= i/- *tp[(£ + A) * qP X [(ft + §) x Q P] + 2 (ft + fi)-^ x[fix(0 X£)] 

+ to * x qP) x C£ x (I x . 0 P) 1 } 


Equation (D.2-56) continued 
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Now 


«“fl*i*n + e*i-o + i0*i-e+(fl + e)*A*{p[j pxecn* p)]} 

L “ - — — ~ — l — — — — — ^/b 0— 0 — 

+ j h *(p[(! -(8 x J) x (ft x (0 x P))]} 

2JS - - 0- - 0- (D. 2-56) 

fi * Cp[qF x 9 (SI- o£) ] > -/ ? *tp[£ *<(£•<>£ - " C&)-*- ~1 2 ( <£" 0 £>i x ^ 3} 

= e_ x i • £ + t tr(i) £ x e 

- e X I -ft + i (I*£ x 0_ + £ x I* e + (I*£) x 9) 

(D. 2-57) 

where TR(I) denotes the trace of the array. Inserting Equation (D.2-57) 
in (D.2-56) and expressing the result In the ji basis yields 

t = ± a m + q m + ~ e ie + j- or+ e xifte + me + eift) 

+ I 7? * { P o P) ~ ^ ® 0 P> (D.2-58) 

Applying Lagrange’s equations 

s} “ [ in + 10 + i I??0 + i Wl 9 -i (in) ~ e ] 

- [-ifiie-yifie+i (i^re - y fi m - y ifi ft + y(ift)~ftj 

- !e l/s * {S[eT o p ^ ^ o p e]} 

* fi iq + 10 + [ia + fit - (ift)~]9 -/i *(cC 0 p ft ft O P0]} 

(D.2-59) 

Now 
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/y <y /v< <-v /v ^ 

Q p n n Q pe - ft q p 2 q p e - <n Q p) 2 .p 9 


rw rv *v> 


= fl q pq 0 p e + (a fi 0 P)'' 0 P9 - ftcfr 0 p)~ 0 p e 
= 0 pn o P 6 + Q 0 p ) o p0 “ n o P 0 P9 + %* Q p 0 


/•o rv /-^ 


- ( 0 pen» 0 P + n^p - + 0 P sT 0 P)e 

- ( 0 pem o p + aa 0 P2r 0 P- 0 p( 0 pnr + 0 ? 0 p $e 


= - ( 0 RnSt 0 t + ff(. [ 0 ? 0 %]~+ 0 P 0 1^e . (D.2-60) 

inserting this in Equation (D.2-59) yields 

Jl M(t)* { o p “ ^ 1Q + ^ + [I ^ + Si “ (inr'je 

+ fs * {p[( 0 ^f^ 0 P + Vo*»~ e - a^pfT 0 ]} (d. 2-61) 

or rearranging 

SH(t) { CT 7 s *<P t ( o p0 ) ^ 0 P]} ■ ^ I« + 16 + [m + fix - (I^p]6 


- Scmr'e + S ine 


(D.2-62) 


Now recognizing that to zero order S = Q Q P this can be written 

^M(t)*^ ( o p “ o pe) ~ s ^ * 0 in + 10 + [in + St - (in)He 


+ [ff 1 3 - S' (in)Te 


Recalling that 


<£ - *V - o p0) " f < 0 p + o P9 > 


(D. 2-63) 


(D. 2-64) 


noting that 


i -y3M(t)* { ^ x $ 
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and recalling Equation (D.2-11), the steady and deviational equations 
become 

T = H I Q 


10 + [l5 + Hl - 0 


+ [ft I ft + K - ft (Ift)~ - 


0 

- T 


3 

- T 


0 0 
0 0 
0 


]0 = 0 


(D.2-65) 


which is identical to the Euler ian formulation in the a basis. 
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